{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas\n",
    "import matplotlib.pyplot as plt\n",
    "from mpl_toolkits.mplot3d import Axes3D"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.seed(1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Multipvariate Normal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(1,) (1, 1)\n",
      "(100,) (100, 1) <class 'numpy.ndarray'>\n",
      "1.0870764500858212 3.1255898724985456\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.PathCollection at 0x7f0e32cbdeb8>"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAD8CAYAAABXe05zAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEc5JREFUeJzt3VuMZNV1xvFvMcPFYGywaSEHmAyRERFBcoNKFhYW7kAcAUZ2HvxglIsTIc2LlUDiyILkIcpDHixFjolkIY3ANk4IuYwhQShBcQgDshRP0oNbeMLgGIgdhoCnkcPFPARwVh7qFNPT09V1qursOnvt/f9Jrenqrqle++xd31l16lSVubsAAHGc1HcBAIDpENwAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDB7Exxo+ecc47v3r07xU0DQJEOHjz4krsvtblukuDevXu3VldXU9w0ABTJzH7Q9rocKgGAYAhuAAiG4AaAYAhuAAhmYnCb2cVmtrbh61Uzu2URxQEATjTxrBJ3/66kZUkysx2Snpd0f+K6AABjTHuo5BpJz7h769NWAHRkZWX4hepNG9yfknTvVr8wsz1mtmpmq+vr6/NXBgDYkrX9zEkzO0XSf0v6OXf/4XbXHQwGzgtwgI6MuuxHHx3++5GPDP/dv7+PapCImR1090Gb607TcV8n6fFJoQ0ASGual7zfqDGHSQAkNOqsR503nXb1WnXcZnaGpI9Kui9tOQCASVp13O7+uqT3Jq4FwHbotNHglZMAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDBENwAEAzBDQDBENwAEEyr4Dazs8xsn5k9ZWaHzexDqQsDAGxtZ8vr3S7pIXf/pJmdIun0hDUBALYxMbjN7N2SrpL065Lk7m9IeiNtWQCAcdocKrlQ0rqkr5jZt83sTjM7I3FdAIAx2gT3TkmXS7rD3S+T9LqkWzdfycz2mNmqma2ur693XCYAYKRNcB+RdMTdDzSX92kY5Mdx973uPnD3wdLSUpc1AgA2mBjc7v6ipOfM7OLmR9dIejJpVQCAsdqeVfKbku5pzih5VtJvpCsJALCdVsHt7muSBolrAQC0wCsnASAYghsAgiG4ASAYghsAgiG4ASAYghsAgiG4ASAYghsAgiG4ASAYghsAgiG4ASAYghsAgiG4ASAYghsAgiG4ASAYghsAgiG4ASAYghsAgiG4ASAYghsAgiG4ASAYghsAgtnZ5kpm9n1Jr0n6iaS33H2QsigAwHitgrvx8+7+UrJKAACtcKgEAIJpG9wu6R/N7KCZ7dnqCma2x8xWzWx1fX29uwoBAMdpG9wfdvfLJV0n6TNmdtXmK7j7XncfuPtgaWmp0yIBAMe0Cm53f77596ik+yV9MGVRAHqwsjL8QvYmBreZnWFmZ46+l/SLkg6lLgwAsLU2Z5WcK+l+Mxtd/y/c/aGkVeFEo05o//4+q0CJRmvr0UePv8xay9bE4Hb3ZyV9YAG1AMB47FDeNs153OgD3RBSG62lWtdWwHET3CMBJw+oAs3LCQju3NXeDWFxaltbgXcIBHfgyZtbTWNFXDQvJyC4o2CxAt0KvEMguANP3sxqfpSBuFifbyO4AdQt4A6B4B4JOHkzq/FRRmpsSywQb+vaN94fAshTxvdNOu7abOwM6Q7nx/MF6AHB3Rfu8Bhh7vMy631zgfNIcNeCHUUaPF+AHhDcfUlxhyc8jomwLdiZ5mna+2YP80hw14LOMC22JxaI4O5bl502nVusbcHONG9t56OHeSS4a1NyOKysSGtr0vJy35UASRHcJaBzO2Z5eTj+SNsiQo2YbIHzSHAjvq0Oj9B5xxVpp9sTgrsPqRYmC/2YUecNFIjgRnwcKsrPLHMR6YnlnhHci8TCBNCB1sFtZjskrUp63t1vSFcS0Jh2xzbvDrDN32Nnu715mhMeObU2Tcd9s6TDkt6VqJbysTCB+DK4/7YKbjM7X9LHJP2RpN9JWhEmy2DhJLXoQ0pt/h6Hudrpojlhm07UtuP+oqTPSTozYS31YGEC8WS0854Y3GZ2g6Sj7n7QzFa2ud4eSXskadeuXZ0ViA0yWjhJLfqQUpu/x2Gu6bB9kmrTcV8p6eNmdr2k0yS9y8z+3N1/ZeOV3H2vpL2SNBgMvPNKgegI/dgy2nlPDG53v03SbZLUdNy/uzm0sSAZLZyFWPT4pjnzAegR53EDqdVyiGuc0sabwTimCm533y9pf5JK0F4GCwdAf+i4oyutmylRbYe4Rmp/pJHQSX0XENbKyrGFWDu2BbBQdNxR0c3E0+XczDPfOZ1mGUVmYyC4pxUxMFPVGHFbAAUguKMqqZtBe/PsLPva0UZem5k2J/UG96wTME1g9j3JqRcdO4/usS3RQr3BXQru4HWZZ2c56/8tcWfSdkyZNif1BXdXXWibTrvvh1eLWnSZLObQclkzCKG+4AZKME+gT9tpl7QzmXVMmY25vuBeRBea28Orvv8+JsttzSBr9QU3pkOQ1KvEnUkhY6o3uBcxYUEXRTLB7ywzWfTnZk6LJypDqje4sb0Sj29iNiXOefAxEdxIr8adQO5jnra+3MdTGYIbWyvkWCBQIoIb6dW4E8h9zNPWl/t4KlNOcHe5oFIszqgLPlq9QAXKCe5UogZujmrchrmPuctXDGNh4gd3l0+abHVba2vS8vLM5c1cHzsM9CnS+otUa0fK/wSceT6dZW1NeuWVYejyKS/IFWuzOvl23H08abLxtkad9qhT7uI2N14eh9Ou0JVZ1k6k9Rep1o5NDG4zO03SY5JOba6/z93/IHVhc5t3UpeXh9etaDEgmIqDK5kg27BNx/2/kq529x+b2cmSvmlm/+Du30pSUQ7v3pXyjafaXi/IAkKG5gn0SOsvUq0dmxjc7u6SftxcPLn58pRFdaKrSa1oMSCYcWu8pOPdiwrlYI9eWh3jNrMdkg5Ker+kL7n7gWQVVbwXPU6t48b8urgP5f5mVxtVeF9pFdzu/hNJy2Z2lqT7zexSdz+08TpmtkfSHknatWtX54XOrMJJRWVKXOOL7oCDNYxTnVXi7i+b2SOSrpV0aNPv9kraK0mDwWD+QymZbzhUJOoraSPch8YFNLbV5qySJUlvNqH9DkkflfT55JUBqQTpqkJItS376oCDrIk2Hff7JN3dHOc+SdJfu/uDacsCMpDi4XqwJ8GmsshPnq9cm7NKnpB02QJqAdLqOjRrDpsUbw+xlRq3bQv5vnKyjZrvOClF2q7RPvS51A5z89tDSLN13mglXnCXtuAjib7tuwrNkg93tJXi7SFSK2ie4gW3NFwoKyt133FSiBRIi6w1ym32hbeHWLg4wb35jrq21lsp1Zn18wlzvQPPW1ephztmEWHskRqSluIE92bLy8ceogWegKxM8xLq0aOevrY94Zkf5mBh4gT3VndUTtZfjLbbfprOpoTAjVx7TQrcyccJ7q0UMAFZ2m675vb8Amugnb7nqVQ9bdd4wc3C68+kRzltOpsCjzciiILWWLzgLlmEECvwYWfR2FGm0fN2LTO4WZz92m67E/zA3MoM7nFyDYuU74mR+m0wkbead5TRXlU7hbKCm4eFcTAndeI+2Ymygnuc3AO9y7137mNFP2qa/+ivqm2hrOCu4TP4gIhoKDpVVnCPE+U4X63vPhepVuQv4n1gSmUGd6QPOgVqUEGYLlKM4O5qsmtaLBHGysNnzKrNWil4HcUI7lxt9ykg0RYNoYlFYH11Iu/gpiMrW6kPn0sbT07IBEm5B3fuNn8KiDTfxzf1gTsCEE7ewV1qR4bjlTKv7ATTIxMk5R7cUWxcPNEWFHcEIJyJwW1mF0j6mqRzJbmkve5+e+rCjkOYIAJ2gotT+bZt03G/Jemz7v64mZ0p6aCZfcPdn0xcW0xRF1TUuoEKTQxud39B0gvN96+Z2WFJ50kiuIGtsBNEYidNc2Uz2y3pMkkHtvjdHjNbNbPV9fX1bqoDAJygdXCb2TslfV3SLe7+6ubfu/tedx+4+2BpaanLGgGUaGWFN4CbUavgNrOTNQzte9z9vrQlZYyFFhPzhsK0OavEJN0l6bC7fyF9SQCKxvnuc2tzVsmVkn5V0nfMrHl5oH7P3f8+XVmZYaHFxLyhUG3OKvmmJFtALQBqwPnuc+OVk22w0GJi3lAoghtAP9iRzozgngYLLSbmDYWZ6gU4AID+EdxASThnvQoENwAEwzFuYCvRzkThnPWq0HGjbhxaQEB03MBGUTtXzlmvCsGNvCwqeMYFNBAAwQ1sFL1zjVYvZkJwIw+LPkQRJaBzrw+9ILiBrRCUyFhewU13Ua9pOuAu10muay3qk6RYCE4HBIBg8ui46S7KNe1ctum0a1gnUY7Boxd03AAQTB4dN91FeVJ0xzWukxrGiKnRcQNAMHl03CN0F+VI2R2zTlA5Om4ACGZix21mX5Z0g6Sj7n5p+pJQFLpjoHNtOu6vSro2cR0AgJYmBre7PybpRwuoBQDQAse4ASCYzoLbzPaY2aqZra6vr3d1swCATToLbnff6+4Ddx8sLS11dbMAgE04VAIAwUwMbjO7V9K/SLrYzI6Y2U3pywIAjDPxPG53v3ERhQDbqun9SYAJOFQCAMHk9V4lwGY1vQc30BIdNwAEQ8eNvNX4HtzABHTcABAMHTdioNMG3kbHDQDBENwAyrGycuz5kIIR3AAQDMe4AcRX2fn+dNwAEAwdN4D4Kjvfn44bAIKh4wZQjsI77RE6bgAIhuAGgGAIbgAIhuAGgGAIbgAIhuAGgGAIbgAIxty9+xs1W5f0gxn/+zmSXuqwnAhqHLNU57hrHLNU57inHfNPu/tSmysmCe55mNmquw/6rmORahyzVOe4axyzVOe4U46ZQyUAEAzBDQDB5Bjce/suoAc1jlmqc9w1jlmqc9zJxpzdMW4AwPZy7LgBANvIJrjN7Foz+66ZPW1mt/ZdTypmdoGZPWJmT5rZv5vZzc3P32Nm3zCz7zX/nt13rV0zsx1m9m0ze7C5fKGZHWjm/K/M7JS+a+yamZ1lZvvM7CkzO2xmHyp9rs3st5u1fcjM7jWz00qcazP7spkdNbNDG3625dza0J8243/CzC6f529nEdxmtkPSlyRdJ+kSSTea2SX9VpXMW5I+6+6XSLpC0measd4q6WF3v0jSw83l0tws6fCGy5+X9Cfu/n5J/yPppl6qSut2SQ+5+89K+oCG4y92rs3sPEm/JWng7pdK2iHpUypzrr8q6dpNPxs3t9dJuqj52iPpjnn+cBbBLemDkp5292fd/Q1JfynpEz3XlIS7v+Dujzffv6bhHfk8Dcd7d3O1uyX9Uj8VpmFm50v6mKQ7m8sm6WpJ+5qrlDjmd0u6StJdkuTub7j7yyp8rjX8gJZ3mNlOSadLekEFzrW7PybpR5t+PG5uPyHpaz70LUlnmdn7Zv3buQT3eZKe23D5SPOzopnZbkmXSTog6Vx3f6H51YuSzu2prFS+KOlzkv6vufxeSS+7+1vN5RLn/EJJ65K+0hwiutPMzlDBc+3uz0v6Y0n/pWFgvyLpoMqf65Fxc9tpxuUS3NUxs3dK+rqkW9z91Y2/8+GpPsWc7mNmN0g66u4H+65lwXZKulzSHe5+maTXtemwSIFzfbaG3eWFkn5K0hk68XBCFVLObS7B/bykCzZcPr/5WZHM7GQNQ/sed7+v+fEPRw+dmn+P9lVfAldK+riZfV/Dw2BXa3js96zm4bRU5pwfkXTE3Q80l/dpGOQlz/UvSPpPd1939zcl3afh/Jc+1yPj5rbTjMsluP9N0kXNM8+naPhkxgM915REc2z3LkmH3f0LG371gKRPN99/WtLfLbq2VNz9Nnc/3913azi3/+zuvyzpEUmfbK5W1Jglyd1flPScmV3c/OgaSU+q4LnW8BDJFWZ2erPWR2Mueq43GDe3D0j6tebskiskvbLhkMr03D2LL0nXS/oPSc9I+v2+60k4zg9r+PDpCUlrzdf1Gh7zfVjS9yT9k6T39F1rovGvSHqw+f5nJP2rpKcl/Y2kU/uuL8F4lyWtNvP9t5LOLn2uJf2hpKckHZL0Z5JOLXGuJd2r4XH8NzV8dHXTuLmVZBqeOfeMpO9oeNbNzH+bV04CQDC5HCoBALREcANAMAQ3AARDcANAMAQ3AARDcANAMAQ3AARDcANAMP8P1KvbS9HeYaMAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "### 1D\n",
    "\n",
    "mean = np.array([3])\n",
    "cov = np.eye(1)\n",
    "size = 100\n",
    "\n",
    "print(mean.shape, cov.shape)\n",
    "\n",
    "X = np.linspace(0, size-1, num=size)\n",
    "Y = np.random.multivariate_normal(mean, cov, size=size)\n",
    "\n",
    "print(X.shape, Y.shape)\n",
    "print(np.var(Y), np.mean(Y))\n",
    "\n",
    "plt.scatter(X, Y, c='r', marker='+')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(2,) (2, 2)\n",
      "(300, 2)\n",
      "[[ 1.97364192 -0.23950527]\n",
      " [ 3.60959539  0.87287921]\n",
      " [ 2.96059751  1.41392441]\n",
      " [ 2.62376799 -0.3208868 ]\n",
      " [ 3.36431546 -0.72461769]]\n",
      "3.330564326051448 0.9197547008358371 0.9589411063033126\n",
      "1.4445657994064762 2.9909216518808908 -0.10179005306793812\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'Y')"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAEKCAYAAAASByJ7AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAGuxJREFUeJzt3X2sJXddx/HPl20F7OKuphuotOuSQAwJMVtyUzUgvfJk1SqikuADRkzYEMXgQ0SxUSSGpIppJILETVssWiEINBAsSIl0cf/g4S5eoLQkKhYoqekS05WFROjy9Y8zw507e2bOPPxmfr+Zeb+Szd57z8N855wzv+/v+Zi7CwCAR8UOAACQBhICAEASCQEAkCEhAAAkkRAAABkSAgBAEgkBAJAhIQAAJJEQAACZS2IH0Mbll1/ux44dix0GAEzKmTNnvuLuRzbdb1IJ4dixY9rZ2YkdBgBMipl9ocn96DICAEgiIQAAMiQEAIAkEgIAIENCAABIIiEAADIkBGButrdX/4CWSAgAAEkTW5gGoEbeKjh1av/vd98dIRhMES0EAIAkWgjAfOQtAVoG6IgWAgBAEi0EYH5oGaAjWggAAEkkBABAJlpCMLPHmNnHzexTZvZZM3ttrFgAAHHHEP5P0rPd/byZXSrptJm9390/GjEmYDjM/kHioiUEd3dJ57NfL83+eax4AGDpos4yMrMDks5IerKkN7n7x9bc54SkE5J09OjRcQMEQmAFMSYi6qCyu19w9+OSrpR0jZk9bc19Trr7lrtvHTmy8TuiAQAdJbEOwd0fNrMPS7pO0j2x4wGCYgUxJiLmLKMjZnY4+/mxkp4n6XOx4gGApYvZQrhC0m3ZOMKjJL3D3d8XMR5gWLQMkLiYs4w+LenqWMcHAOzHSmUAgCQSAgAgQ0IAAEgiIQAAMiQEAIAkEsL0bG/vLXACgIBICADCo+IySUlsXYEG2CANwMBICADCoeIyaSSEqWCDNGAeEr6GSQgAwqHiMmkkhKnhAgOmaQLdaSQEAOElVMihORICAIxhAt1prEMAAEiihQAA40qwZZCjhQAAkERCAABkSAgAAEkkBABAhoQwRewkCWAAJAQAgCSmnU5LKkvfE15YA6A7WggAAEkRWwhmdpWkt0p6vCSXdNLd3xArnkmIvfQ9lRYKgEHE7DJ6RNLvuvsnzexxks6Y2V3ufm/EmABgsaIlBHd/UNKD2c9fNbP7JD1REglhk1g18tgtFACDSmIMwcyOSbpa0sfW3HbCzHbMbOfs2bNjhwYAi2HuHjcAs4OSTkl6nbu/u+6+W1tbvrOzM05gQF+0pJAIMzvj7lub7he1hWBml0p6l6TbNyUDAMCwYs4yMkm3SLrP3W+KFQcQHLOxMFExWwjPkPQSSc82s93s309EjAcAFi3mLKPTkizW8YHBMBsLE5XELCPg29i4D4iGvYwwfanWxFOLB9iAhIA0MBALREdCwHSlnERSimUqeM2iIyEgDQzEAtGREDBdKSaRlFstqWr7mvGaDoaEgIvFvOBSv8gpjDBjJARMX+jCuU+hn2KrJXVNXzNaX4MjIWAPF1w1XhssAAkByIUs9EkU7W16zWh9DY6EgD1LvOCanusSXxssDgkByFHoTwPvy2BICLjYEi64rt1Dc35tSISLR0LA9JU3w+tboFEgYqFICFim0N1DU65d9x1Mn/K5Yx8SAqarXJCV/04BBbRCQsCyhWoZTHl9QtfW0hzOHfuQEDBd5YKs/PcmKMSAbyMhIH0pF9pzmqraNvY5nTskkRDC4aKIp8uumKG7O3Z3uz0OSAgJAemaUh/18eOxI4gnxfcDnZAQ+ppSobUUTd6TNt0dY7Y0gIhICEhXjD5qCnQsWNSEYGa3Srpe0kPu/rSYsXTGwFp62rwn627b3d0/cylUSyO0FD9zKcaExmK3EP5W0hslvTVyHFi6ctdPPki85LEBLE7UhODuHzGzYzFj6KSupojuQtcu2z5PcabQuXPSoUOrn6+9dvPzdYm56/mmOG6RYkxoLXYLAVM35IU/ZiFz992r59/dXSUDidYBFif5hGBmJySdkKSjR4/GDaZcQB0+vCo0qAX1U1Xwj62YFIZ6X4vP3zXRpbgxH2Nps5B8QnD3k5JOStLW1pZHDge5MWrvMQqZPCmMLR/IpiBFRMknhKTkhcWhQ6tuhbxrYSoXcqq1t9i1y/Jx61Y+d209rNuZ9dChvefqkoRS3Jgvtc8WWok97fRtkrYlXW5mD0h6jbvfEjMmNDRmIT7kcxcL+eLfhj5uLm8ZMBiLBMSeZfQLMY/fSezabBdTmQESo2WQDyKfOlVdSy8PNp861X78aN3OrLEHraf4WcagltdlNLcP/9CDoJvEeh2L72OX97RcyEvS6dPSwYP7uwJDxVilagvvIY85hWMgiuUlhFCmdDFQE1wpn39xpk++5qD4t1zx9QuVfGPPqCpa6ucBF1lOQphKt0lTIboxpmjd1N9yjb7J+dcV8n0/G30+a20fU04oQ36+53YN4SLLSQjov5p2qvv1bCrImibREOddlXSGbClQcKOh5SSElLtNusQ0RDdGDF0XZR0+vPr/4YfXP0fT5113e99pp10+a20fU05y+fYaTbbZ6CrlawhBLCchoJ1QXTOhjr/umOfPN6tZD1UL7/M69H3tmnxDG108aGl5CSGliyHEBdt2Zk3bxwyl6zYg+ePy5FTV/x+yEKwaryne3qalUafpY8pTVsd8T1P4/GAQy0sIaGZd98DYW0hIewVvuYssd+7c/m0fmtScQySL8+elCxdWP4+57UQ54eWzo9ahiwctkRBiGuuCDV1rDrXFQd02IOWWQDkxVe1KOkTNufg+nT69+rm4XqGclMYsgGMvbpsiEmQlEgLqFS+atn3xQ059PH58byFZXggXt4BYV3Puk4DzJJS3CsrGKpj7DFh3QeG5KMtJCCl/sIeOKYX59evUxbVpeubdd+/vxy/rW0Cvi6m4aC1PSJL0zGeu/m+7L1HMz2TK18NQGGTfaDkJAeMZamC3+JxVs57abBfRRrnbqNhaaDJuEdrQhRiF5yLNPyHwwd7T95yHGvMYc6+forr1C3WflwsXVoPKRcU9kYrbWtcdt5jcpNWaiqENmaxTv6YYZN9o/gkB4xviwms666nvauwm95NWSeHQoVViWLf/0dRReC7S/BMCH+zwUn8Nm/bdrxuArvu8lLuG8lZC25lF+X0uyS6/MRf8FWd3bVr3UVyZXXW7NL3Wd+rxRVSZEMzsTkm/7u73jxcOZmWI2UVVrYG2hVF5LUOX5zh4sNn9yvLjVM1YSsmQW6JMJYEsSF0L4S2SPmhmt0n6c3f/5kgxDYMPXRgpXcTlWJrO8slryeXvPyjfpyyvKRdnGtXdv6m8hTL2QHH+t/Jxm9b8aX3PTmVCcPd/NLP3S/ojSTtm9neSvlW4/aYR4kNosVYed7WpcNrd3T/jaNOMn3WFYpPuk+Ixy6un2xpyoHwKptrVtACbxhC+Ielrkh4t6XEqJASUzP1DndJFXG4JHDiwd9vx4+12f83v30boRWhjvYZNa/Rta/5z/cwvUN0YwnWSbpL0XklPd/evjxYVwgu9e2lVIRp6MVbd7KJc3hdfHBxeF1+o1tEY023njK6mZNW1EG6Q9CJ3/+xYwUxSSjXnIY2xbUWbGun29qplsG4MYFMNftPsmbbx9H1MH12O1zYZYzHqxhB+ZMxAMLCm8/g3yWveVV9EX7WXUN/Eue5+Bw+uL9Srvq+4HNtUCry5VzKQjPmvQxhal+bvlC7wcqzlmTZVqhJGrtil0yZZrBsUblLbXxdb1QrhLslrXZdcqFlITY43pc8UkkVCWJq6Vb1t9gGqum9eyJa7cULNrFk3NlAeQG4bW6ryrq0xF65h0aImhGzg+g2SDki62d1vjBlPL21aBlOo1XWNddP3EVR1OTXtwurzfQf5Y+tiKP7fdvVxvgI4/36HYmEeer+iqi7AdesKgIaiJQQzOyDpTZKeJ+kBSZ8ws/e6+72xYlqsdYX/pm6YplMWy8codzkVa/x10z/rtpuoi61YYJbXEaQsjztPKDHHPVKuuEjpxzchMVsI10j6D3f/vCSZ2dslvUDSfBPClKbb5d0wocdG6rqO2owD9FkLkNfSm87Hb2rd+1vuphpip9i8ZTCFlieSFjMhPFHSlwq/PyDpByPFMq6h988PMde/y/HqbisXVuVum03rItrGOIXuubqYUqg8pP4aph7fBCU/qGxmJySdkKSjR49GjiaQKX0PbpOLq65wr+raaTpbaUhDFRzF523aGglxvKVuhYFgYiaEL0u6qvD7ldnf9nH3k5JOStLW1paPE9pAhq7RDDHXv8nx6mbtbBoEXte90qTLaZMxathdV2W36UZqG3eenENUOlJopdRJPb4JipkQPiHpKWb2JK0SwYsl/WLEeNBX3hp4+OGLL9K62mu+N1HZnC70Mc/h3LlVpSDW6zen963ODM8zWkJw90fM7BWS/lmraae3zn6bjKFrNF2fv2s8bboq6gZuDx9eP9++y1hL1RTSkKoWoZVfj6oV01Uzh/p8LsprFqT2C/aqpF7gpR7fhEQdQ3D3OyXdGTOGVmZYIwiqbtHbOlVdTpvWCSBNSxnknfF5Jj+oPEtDf3DatgyafrCrbg91PsVtrIuatBTGvEjXLULLj7luvcOmFeEhYs9jKrYS8q29QyxWm1Ghh2okhCZmXCOIqlxwVg1Od+32aFsYtn1f8+9UrnqesRWTQnkdSVtdk1LXx07JjM+ThLBkTT/YY++pU95iYt1xqsYKirXyqsIwRPzrCt9c060xis81RFwhFqtRGVoUEkITM64RJKFqrn5dgV4eMC0nrTZfANS166y4/UYx1tCb53WZPty3ZdBki5C64xefqxj3nK6hOZxDCQlhyZpcnDH31KlqGeQFf3FqZRObZv90VdWl1fa1CvmahqrE9NlMEJNDQmhjiRfD0NtsFDXZjmLd1EppVXAVk9a67yKoSgBtCs/ifcrPV3WcKXXXVL0WfRNv3qKj6ylpJIQlalPQhKohVvW3t1UsVMoF/yZtxkzqbu96vLEU4w9V4FJwxzNi8iQhYL9y4dWnL3ksobpH6grww4dXs4ouXNi7T1WtNxdqK5G6xw29yLHPY+c+hjBDJIQl6lKAdp36We7mOXXq4pW9XXWNKdTgcmqmHj/2i/B+khCwUv7wXXvt/v/HLFQ2LYAbeoFc0e7u/hlL0t4CuvLsqFzo17CuZdDmS43Gti5uklPSSAhLNuZMoVBjCENal3DyrqIp6PKlRkhXhOnuJASsDPnh296WTp9udj9pfa237T5JoRRr3FUJbdPWFENgbQwGQELAOA4eHK51ELpQzGvZxZW+hw4N11IYIn7Mx4jvJwkB+4VuGbQZUC5uP5EXvn339g9V2OYthSbPM2aBTOGPgEgIc9Z2odWY6qYk5glk06K4IWdhjNElw6wgJIaEgOHUDSg3WayVb99c19VUlTRiFrYU7JgoEsIcNSkMYxWYdcddVyuvSxx1O46GMvTA8Pb23orrPsciCSEAEsJcpFwgrEtEXR5bVLXjaH7/GLNwYrZKxtxzCrNFQpijJoVhrGmLbWJrI6UFWU2Uk0f+t7bbU5THXVKuGCB5JISpa1MrHaqw6Pu8XR7fNKH1Pdc2sfVJsl1XGZdbBiFbCiSXxSEhzFlqUyS7HreqYCpv1zClAqzJCu4myT4/9/w+U2spISkkhKlrsyNm6L7tvs9bfnz+9Zf5HkFNtBmMbRNfeQ1FHluT47WJZ92X/XRpKRW/OrQvpsMuFglhiaYyAFlVMOXKiWRTP3oK512Orbh9dlmbLqgxWwYkiNkiIcxF3cXZZu5/l2N2LSDKtdtygR5K24HX8v2LP/ddNV3s4go5sF83k6vr+0LBvzhREoKZvUjSn0h6qqRr3H0nRhyz0fTC3d29+Etemjyuyf3K0z5D2FQwlf/etGUQo6XQp3solQKZrqTZi9VCuEfSz0r6m0jHX6Z85W+x5hvCpgVkm5S/VyB0AdN24HVTi6pry6D8XdBt9kdqehzp4i6ovi04LEaUhODu90mSmcU4/Hx0mXJa/CL6qvu2ef4xao1Np5Ruul/IgVep3bmu+y7oqRW4dCXNXvJjCGZ2QtIJSTp69GjkaNBYl+6ooQuYcsugzxqG8nhAnSZTTLsox1/VsqHgRkODJQQz+5CkJ6y56QZ3f0/T53H3k5JOStLW1pYHCm8e2tTYutTu+qx4Dj0wHELIlkHX8YBUWwZTHNNAcIMlBHd/7lDPjYQ17UJKYTfSdcduMgOpbjygzlBjBVUtBaCl5LuM0ECXhUyhnz/llkFoocYDUuiLHyMxT+U5EW3a6Qsl/ZWkI5L+ycx23f3HYsQSFB/S9nsMxUgg62Lc3l7921QwFgeoz5+P0wXE4C4GEmuW0R2S7ohxbKCXdYvW2q6/6FMrD50EhkwuQ7Q+WAsxKLqMQkjtQxrywuv6HG0f13TMYajprFOsdYecqQSIhADsLxQ3DRCHSBxdnmPoSscQCbA83TbEcaaYuCeEhBBCKh/SEIXG2K2dpttT9I0n9nuTmtRatUgCCQHLVbXVw+nTq0Hjum24QxScXWaHNd2YL5WCvTjm0mdjwLJUzm9mSAghxf6QblokFmox2jpDjTn0bX11qQlfuLCaQdRksDi1AripVFq1SAoJActV7uMu7gR74cLmHVzbbF9Rflzx+G1sahmk1gVE4pkUEsIclS/CLoVE0y6JsQqivi2PTXHlO8EWk0JdQV+1fUXfeMc2lTgxChICUCwUm2xA13X7iiGTZ+o18dTiwVokhDkLWUhsKsxSLYhyTeMq7mdUp7x9RS61Lpsl4rXvjISA6YrVPZUnjTwRNB1DGCN5LqkQpOAPjoSwBGN0SUz1ouy6ViNfaFVMBk12S8VwUh1YnxASAqYn9tbZXbe/zlFA9VP1/qM3EgLaGaowG7s21yWprEsGp0+vZiWFXHSFbqYynpUwEgLiazufP+aFXxxMPnBAOnhwf4LA8Cj4B0NCQFxV8/mHvsj7ruIuTk2lYEoL70NnJATE07c/fsgLv8lCNgqeuHj9gyMhIK5QX0fZVYjtIyiYMBMkBMRT1QUT05SmLqYcGybpUbEDACSlkQyAhaOFgPhSSgRDzWAJ+XxTasWkjNftIrQQlmR7e9mLeIrnP/fXYu7nh0HQQgDWCd0yCFmbZx5+P7SwKpEQlmDpF0D5/A8f3pvqOrfXYunvNXohIQBDGrI2TyHfDS2sSlESgpm9XtJPSfqGpP+U9FJ3r/lGc/Sy9Atg3fnP9bVY+nuNXmK1EO6S9Gp3f8TM/kzSqyX9fqRYgOFRMKeH9+Qi5u5xAzB7oaSfd/df2nTfra0t39nZGSEqAJgPMzvj7lub7pfCtNNfk/T+qhvN7ISZ7ZjZztmzZ0cMCwCWZbAuIzP7kKQnrLnpBnd/T3afGyQ9Iun2qudx95OSTkqrFsIAoQJoizGKWRosIbj7c+tuN7NflXS9pOd47H4rYEoojDGQWLOMrpP0KknXuvvXY8QAoAPWOcxarFlGb5T0aEl3mZkkfdTdXx4pFmAaKIwxsCgJwd2fHOO4AHpincOssVIZmAoKYwyMhIBxUZjNA+/fLJEQgKmhMMZASAgYBwOiQPJSWKkMAEgALQSMgwFRIHm0EAAAkmghYGy0DIBk0UIAAEgiIQAAMiQEAIAkEgIAIENCwGoqaD4dFMBikRAAAJKYdrpsbCfRH68ZZoQWAgBAEi2EZWM7ie5oXWGGaCEAACTRQoBErbYLWleYIVoIAABJtBCAfmgZYEZoIQAAJJEQAACZKAnBzP7UzD5tZrtm9kEz+94YcQAA9sRqIbze3X/A3Y9Lep+kP44UBwAgEyUhuPv/Fn69TJLHiAMAsCfaLCMze52kX5F0TtKPxooDALAyWAvBzD5kZves+fcCSXL3G9z9Kkm3S3pFzfOcMLMdM9s5e/bsUOECwOKZe9zeGjM7KulOd39ag/uelfSF4aMK5nJJX4kdxAA4r+mY4zlJnFdb3+fuRzbdKUqXkZk9xd3/Pfv1BZI+1+RxTU4oJWa24+5bseMIjfOajjmek8R5DSXWGMKNZvb9kr6lVY3/5ZHiAABkoiQEd/+5GMcFAFRjpfKwTsYOYCCc13TM8ZwkzmsQ0QeVAQBpoIUAAJBEQhiEmd1qZg+Z2T2xYwnJzK4ysw+b2b1m9lkze2XsmPoys8eY2cfN7FPZOb02dkwhmdkBM/s3M3tf7FhCMbP7zewz2V5oO7HjCcHMDpvZO83sc2Z2n5n9cJQ46DIKz8yeJem8pLc2WV8xFWZ2haQr3P2TZvY4SWck/Yy73xs5tM7MzCRd5u7nzexSSaclvdLdPxo5tCDM7HckbUn6Lne/PnY8IZjZ/ZK23H026xDM7DZJ/+ruN5vZd0j6Tnd/eOw4aCEMwN0/Iul/YscRmrs/6O6fzH7+qqT7JD0xblT9+Mr57NdLs3+zqCWZ2ZWSflLSzbFjQTUzOyTpWZJukSR3/0aMZCCRENCRmR2TdLWkj8WNpL+sW2VX0kOS7nL3yZ9T5i8lvUqr9T5z4pI+aGZnzOxE7GACeJKks5LeknXv3Wxml8UIhISA1szsoKR3Sfqt0s61k+TuF7Kt2K+UdI2ZTb6bz8yul/SQu5+JHcsAnunuT5f045J+I+uinbJLJD1d0pvd/WpJX5P0BzECISGglayf/V2Sbnf3d8eOJ6Ssmf5hSdfFjiWAZ0j66ay//e2Snm1mfx83pDDc/cvZ/w9JukPSNXEj6u0BSQ8UWqbv1CpBjI6EgMayAdhbJN3n7jfFjicEMztiZoeznx8r6XlquLdWytz91e5+pbsfk/RiSf/i7r8cOazezOyybEKDsm6V50ua9Gw+d/9vSV/KtvORpOdIijJRI9r3IcyZmb1N0raky83sAUmvcfdb4kYVxDMkvUTSZ7I+d0n6Q3e/M2JMfV0h6TYzO6BVBekd7j6bKZoz9HhJd6zqJrpE0j+4+wfihhTEb0q6PZth9HlJL40RBNNOAQCS6DICAGRICAAASSQEAECGhAAAkERCAABkSAhAR9nur/9lZt+T/f7d2e/H4kYGdENCADpy9y9JerOkG7M/3SjppLvfHy0ooAfWIQA9ZFt5nJF0q6SXSTru7t+MGxXQDSuVgR7c/Ztm9nuSPiDp+SQDTBldRkB/Py7pQUmT3yUVy0ZCAHows+NabYj3Q5J+O/tWOWCSSAhAR9nur2/W6nshvijp9ZL+Im5UQHckBKC7l0n6orvflf3+15KeambXRowJ6IxZRgAASbQQAAAZEgIAQBIJAQCQISEAACSREAAAGRICAEASCQEAkCEhAAAkSf8PmnAjytgyUjYAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "### 2D\n",
    "\n",
    "mean = np.array([3, 0])\n",
    "cov = np.eye(2) # 协方差Cov(X, Y) = 0, D(X) = D(Y) = 1\n",
    "size = 300\n",
    "\n",
    "print(mean.shape, cov.shape)\n",
    "\n",
    "XY = np.random.multivariate_normal(mean, cov, size=size)\n",
    "\n",
    "print(XY.shape)\n",
    "print(XY[:5])\n",
    "print(np.var(XY), np.var(XY[:,0]), np.var(XY[:,1]))\n",
    "print(np.mean(XY), np.mean(XY[:,0]), np.mean(XY[:,1]))\n",
    "plt.scatter(XY[:,0], XY[:,1], c='r', marker='+')\n",
    "plt.xlabel('X')\n",
    "plt.ylabel('Y')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$\n",
    "X = [\\mathbf{x_1}, \\mathbf{x_2}, \\cdots, \\mathbf{x_n}]\n",
    "$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Scatter Matrix"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "S = \\sum_{j=1}^{n}(\\mathbf{x_j} - \\overline{\\mathbf{x}})(\\mathbf{x_j} - \\overline{\\mathbf{x}})^T\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 142,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(3, 500)\n",
      "0.03689421074351941 -0.02067275040562437 -0.08732962492218015\n",
      "Scatter Matrix:\n",
      " [[474.83941486  -5.49289844 -26.38804904]\n",
      " [ -5.49289844 514.38666859 -39.1145724 ]\n",
      " [-26.38804904 -39.1145724  496.45740348]]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<mpl_toolkits.mplot3d.art3d.Path3DCollection at 0x7f0e2995b3c8>"
      ]
     },
     "execution_count": 142,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAV0AAADnCAYAAAC9roUQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzsvXmcHAWZ//+uo6+5J5mZJHPkmgwhISQhgSSgsiKCgorroq7HLkRF3AU0IBsN/FAXxQCyEFxRXA7BRSB+FRUV5VgQBSEkXCFAjsnc990903fX8fujuorqnu6Z7p6DDPTn9corPT1TR1dVf+qp5/k8n0fQdZ088sgjjzxmB+LbvQN55JFHHu8m5Ek3jzzyyGMWkSfdPPLII49ZRJ5088gjjzxmEXnSzSOPPPKYRciT/D4vbcgjjzzyyB5Cul/kI9088sgjj1lEnnTzyCOPPGYRedLNI4888phF5Ek3jzzyyGMWkSfdPPLII49ZRJ5088gjjzxmEXnSzSOPPPKYReRJN4888shjFpEn3TzyyCOPWUSedPPII488ZhF50s0jjzzymEXkSTePPPLIYxaRJ9088sgjj1nEZC5jeeSRFrquo2kakUgERVGQZRlRFJEkCVEUEUURQUhrtpRHHu9KCJMMpsxbO+YxDrquo6oqiqIkvDZ/ZxKtpmmMjo5SWVmJJEl5Ms7j3YS0F3g+0s0jYySTrSAIiKKIpmnWazs0TaO1tZWysjKi0WgC0ZoRcZ6M83i3IU+6eUwKXddRFAVVVRPIdjKYfyNJ0rj1ASiKQiwWG7dMnozzeCcjT7p5pIVJtmbqIJlsNU2ju7ubtrY2NE3D6XRSWFho/fN4PKRKX5kEmkykyWRsT1UIgoAkSVbe2CTnPBnnMdeQz+nmMQ6apiXkaQVBSCA3RVHo7Oykq6uLqqoqamtrEQQBRVEIBAIEg0ECgQB+v5/R0VFKSkoSyLiwsBCn05nx/pjXqK7r1j9BEOjt7aWiogK3250yMs4Tch5vI9JefHnSzcOCpmlWGgHGk20sFqO9vZ2enh5qamqoq6tDlmVUVSUWi41LOei6zosvvsj69esJBAIJ/2KxGLIsU1hYSEFBQQIZZ0qWr732Gg0NDbhcrnERtRmVy7KcJ+M83g7kC2l5pIYZOcZiMbxeLx0dHaxZsyaBmCKRCK2trQwODlJXV8epp546Lk+bCuY6HA4HZWVllJWVJfzejIwDgQBDQ0O0t7cTjUaRJMkiYZOQXS5XWrJMRfaAJWdL3iczNWFPVeTJOI/ZQp5036UwNbaKoqBpGmAUvOx51FAoREtLC16vl6VLl9LQ0JBRAS1TyLJMaWkppaWlCe8rimKlKEZGRujs7CQSiSCKYkKKQtO0rHLG9s+tqirRaDThd/YCnhkd58k4j+lGnnTfZUgmW5NU7PIvv99PS0sLgUCAZcuWsWrVqlklHlmWKSkpoaSkJOF9VVUtMvZ6vfj9fg4cOIAkSQkpisLCQtxud8p9TkeidjI2ibylpYXly5enVFPkFRV55Io86b5LkEpjm0xAfr+fkZERDh48yLJly5g/f/4xRSySJFFcXExxcTEAwWCQFStW4HQ6LTL2+Xz09PQQCoUQRRGPxzNOUZEJGWuahtfrtW5Edrmciby8LY9ckCfddzjSNTTYMTw8THNzM7qu4/F4OOWUU96mvc0NyWRsQtM0i4zHxsbo7e0lHA4DUFBQQEFBAUVFRRQUFODxeMYdl1Q3JhP2XHi+8SOPbJAn3XcoJmto0HWdwcFBmpubcbvdrFy5EpfLxf79+9/GvZ5eiKJIUVERRUVFCe9rmkYoFLJkbX19fYRCIQArMvZ4PGiahqZpKfPYE5Ex5Bs/8kiPPOm+w2CSbX9/P6FQiNra2nFk29fXR0tLC8XFxaxZs4bCwkIAqyEhW8w10rAX5OwwyTgYDDI2NkYkEuGll14CwO12J6QpCgoK0pKx/X8TqRo/jh49ahUnTXlbvvHjnY886b5DkNzQoGka4XA4wXymu7ub9vZ2ysvLOemkk3C73QnrMPOX71bYybi8vByv18uGDRvQdd2KjE15WzAYRNO0lGScSk6XTMaqqhIIBCziTk4B6bo+YWScJ+S5izzpznGka2iQJMki4s7OTjo7O6mqquLkk09O2w1mftnzIOE4CIJg5YArKysT/iYcDltkPDw8nEDGyYoKOxnbi5nmNlJt334jte9PvvFj7iJPunMQ9iKOGZmmqr6PjIzw/PPPU11dzebNm5HliU93nnQTMRmBCYKAx+PB4/FQUVFhva/rOpFIxCLjrq4ugsEgqqricrksSZtJqKnOy2RknG/8mLvIk+4cQqqGhuQvVDQapbW1ld7eXpxOZ8bdY+a65hpm6iaRLA/LBoIg4Ha7cbvdzJ8/P2Gd0WjUkrZFo1H279+PqqrjzIIKCgpwOBwp123/P3mf7Y0fzc3N1NXV4XA48o0fxxDypDsHMFFDg4lQKERraysjIyMsWbKEdevW0d7enjHhzkXMNGFM9/oFQcDlcuFyuSgoKGBsbIy1a9cmkHEwGKS3t5dAIICiKDgcjnFmQanI2Fx/su7aJNhkrbH5f77xY/aRJ91jGKbGtru7m3nz5iHL8rgvViAQoLm5mUAgwNKlSzn++OMRBIFAIGDleWdjP9NJq+YqZjrNYj9edjKeN29ewt+ZZBwIBOjr6xtHxslmQam2kWkXXl5rPDvIk+4xiOSGhs7OTkpLSxMinNHRUZqbm4lGoyxfvnxc95goijNOHKb8zGysMJsr7FHZbOWIHXfcAUDs4ounZX1TSS9M5/qdTidOp5Py8vKE92OxmEXGAwMDtLa2Jji3FRYWEo1GLfe3bFui840fM4c86R5DSNfQYJdyjYyM0NzcDMDy5cvHfRlNzKT8S9d1enp6aG1tpby8nPXr11tpjHA4jN/vJxAIMDg4SDAYZN++fSnbcacrMnbccQdiS4v1erqIdyYxVVJP59wWi8WsLjxFUTh06NA457bJbDQzbfxQFIX29naWL1+eJ+MskCfdYwCTTWgQRZHBwUGrONbQ0DDODCYZM0G6mqbR09NDW1sb8+bNY+PGjbhcLjRNIxqNWl4HHo/Hklbt27ePDRs2EIlELDI2GzeAcbKqdN4Is42ZjnRnKh3jcDgs57auri7Wr18PvOXc5vf7GR4epqOjg0gkktIsKJ2NZnIRT1EUy+PC/DlVF146edu7FXnSfRsx2YQG8/F9aGgIVVU54YQTxnVRpYMoitOW09U0ja6uLtrb26msrEyp9U33JTI/k6lzTV6v2XRg90Yw/z4T1zAwUgrv1PTCdGEi5zYzTTEyMkJXVxfhcBhRFCc9B6qqJhDoZFrjVI0fdmnbu0VRkSfdtwGTTWiwR5Tl5eVUVlZSU1OTMeHC9ES6qqrS1dVFR0cHVVVVnHLKKSkbK3L9ktg7wKqqqqz3TaMav9+Pz+eju7vbIgI7CdhtGOdCSsEOU4XydkOSpEltNJPPgXkDNesG6W4g+caP1MiT7iwhk4YGe/eYPaI8fPhw1lHrVC5Q+wy0hQsXsmnTprQypUz2I9tiWjqjGpMITAvKsbExDhw4kFJWlc3Yn1SYjUj3WFZ7pHNuU1XVejoZGhrC7/fz4osvAuNTRW63Oyd/CnvjRywWo6+vj7q6undM40eedGcYdo1tR0cHgDXI0YRZkOju7qa6unocyZktvTMNRVGIRqPs2bOHmpqajLrYJsN0drklE0E4HGb58uWWn67f708Y+2Ov5Gc7EHM2lB8zTeozAUmSrBuiLMu4XC7q6+utpxPzPKRybpusiJqKjCORCGNjYwiCkNHEDzOomazm8XYiT7ozhFQNDaIoEovFrIsqGo3S1tZGf38/tbW1abvHJEmaUc2toii0tbXR29uLIAiccsopuFyuaVn3TLcW67qeNl9pVvL9fv84WVVRUdGkDQczXUib6zljVVWt6zUTG01T3hYMBoHMnNvMNulMtMYAjz32GPv372fnzp0z8ZGnBXnSnWZMNKFBkiTC4TDhcJiWlhZGRkZYvHgxp5566oSPmjMl/4rFYrS1tdHX10dtbS1btmzh5ZdfntYv60yS7mT7aa/k2xGLxSwlRXLDgUnGs3GzmMn0wmw0q9hJNx0mstG0mwWZ8kJT622mKsxiXTokk7HP5xt3vo815El3mpDJhIZYLEZvby/9/f0sW7bM6h6bDNMd6Zr+DAMDA+NIf7oJ/ljMtzkcDsrLy8dpnM3uL7/fj9frxefzsXfvXsukJp1jWC6Y6Uj3WCHddLAX5CZybhsaGiIUCuH1ejOy0fT5fOO0y8ca8qQ7RUw2oQFgbGyM5uZm/H4/Ho+Hk046KasvnJmWmCrMUepDQ0NpI+yZiPDminOZvfurtLSU7u5uVq5cmUDGdsewTL10U2GmI93JIsTp2kauBdZ0SHZuczqdKIpCbW0t4XA4YUp0IBBA0zTrpvj888/T0tLCli1bpmVfwuEwp59+OpFIBEVR+OQnP8m111475fXmSTdHTNbQAOD1emlqagJg2bJlSJJER0dH1hHOVKJP02bQTGcsnWSU+kxEunOFdO2wp4ZS+SLY7RvNhgPTSze5cJQqVzkbzRczbXY0lUg3UyiKYk3RMMk42bnNPA+9vb288sor/O1vf+OWW26hpqaGhx9+OOdisMvl4qmnnqKoqIhYLMZ73/tezjnnnCmTep50s4SmaYyNjaEoCgUFBSkbGoaGhmhpacHhcCR0j/n9/pzSBLmmF3Rd580332R0dDTjdMZ0ezbMVdKdDBPZN5qPx36/38pVQmIVPxQKTVuxMhWO9fRCplAUZcLjZD8P27dvp62tjUsvvZRTTjmFvr6+KalvBEGwCoOxWCyhCD4V5Ek3Q9gbGgYGBlBVlaVLl1q/13Wd/v5+WlpaKCwsZNWqVeMqubmSZ7bdZcFgkObmZoLBIMuWLWP16tUZXyymNGe6MFdJN9dINPnx2L4++zDMoaEhFEWhp6cHj8eToKaYjlbo2SDd2YqmsyHO0dFRSktLEQSBhQsXTsv2N27cyNGjR7n00kvZvHnzlNeZJ90JkK6hweFwWFpBe/dYWVkZ69atw+PxpFxfrqSbqU7XtHkMBoMsX76cYDBIRUVF1vnjfKQ7/Y//9lboyspKJEnC6XRSVVWVshUasEbEZ9IKnYzZIF3z0X+mt5EN6fp8vrQmULlAkiReffVVvF4vn/jEJ3j99ddZs2bNlNaZJ90UmGxCgyzLRKNR2tvb6ejooLKy0jJ/mQhTId2JlvP7/TQ1NRGJRBJsHtvb27OOWnOJdM2uIVPmk/xFnIukO9MwSTGXVuhkMk5lUDNbke6xRuwzpV4oKyvjjDPO4NFHH82T7nTClH2pqpp2QoM53rynp4elS5dm1SKba5Eq3XJjY2M0NTURi8Wor68fZ4Cdy/ayWSYWi9Ha2kp/fz/z589PqCib+ctIJEIoFEo7svxYxdtteJNpK3RnZ2eCW5hJxpFIZFabI2ZyG9lEutOpqBgYGLAsNEOhEE888QTf/OY3p7zePOkycUODCXv3WEVFBfPnz6e+vj6r7eT6JUiOdH0+H01NTWiaRn19/bR66maSDrCT7eLFi9myZQuKolifzywmmYWkrq4uWltbAeORuf7uu3E4HMSuvTatjWCmmMkZaTOJXHW66TwRTOtGU9s6PDxMLBZjZGRkXPddpq3Qk2G2CmmZku50n7Oenh4uvPBCKwj79Kc/zUc/+tEpr/ddTbqZNDSEw2FaW1sZHh62tK3hcJjDhw/P2n6a5GmXoNXX10/6GJWLZ8NERG3vYLPrfJMvdnsxaWhoiOrqakpKSozCy7e+Bc3NqIpC+Moree3CCy2DbTs5ZEIMMx3JzSXDm+RW6O7ubjRNY+HChWknTGTSCj0RZkshkek2zOtwus7b2rVreeWVV6ZlXXa8K0nX1Ni2tbVRVVWF0+kcd2KDwSAtLS2Mjo6ydOlSVq5cmZDTTbakm0mMjo7i8/loaWnJyMDcxHSlF+xkW1dXl7apIhXskbMoijicTiSXC1wu3FVVlG3ahKIoKYnB3pZr/j9bgzbf7vTCVGESoizLaVuh081eSybjiSLN2eg4zHQb4XB4nGfzsYh3FekmNzQMDw9TXl6eUAAzu8dMB6tUcquZNqAx93V4eJimpiacTicul4uTTjopq3XkUhSzk6SiKLS2tk5IttmsDyD6ne/gjHf1RL/zHYC0xGDvBOvu7raGbZqdYKFQiGAwOOfyxTDzUaKmaROSZbpxP/ZBmL29vZa23BwRbyfkYw1er/eYdhcz8a4g3XQTGuwRq9frpbm5GU3TrNlj6e6wU+naMokw3RdO13UGBwdpbm7G4/GwevVqioqKeO6557LeVq6RbiwWo6mpid7e3pzJ1kSqHLFJtpMh1VBGe/PB0NCQJdeDt/xcTWKwS6ycN9xgbHvHjoy2/U6JdLNFumNuJ2OzFToQCLB///6EY55NK/RkyPYYzQXfBXiHk+5kExokSWJ4eJjGxkZkWaa+vj4jh6KpfFnMKDlVW+jAwADNzc0UFhayZs2aKUcT2eZ0FUVhaGiIsbExli9fnhXZzhSJOO65B4DYF74AJOaLzUJeUVGRZSHo9/sZHR21JFaSJLHigQco7uxEliSk738f9f/7/zL6PDOJY5V0UyFVK7Su6+zbt4+VK1daZGxvhTafRuxknO3+ZPsZ5oLDGLwDSTeTCQ0mwfX29uLxeDjhhBPGSXNmCibpmkULcw5aS0sLJSUlrF27dtryUplGunY/3cLCQpYuXcrixYunZR+m0hzhuOcehLjqwXHPPRbxpkI6C0FVVRH/+EcQBCLRKKP9/TTt3WtNm7A/Lic/jqcjxeQbQS6YjfTCTK9fkqSsWqFN20b7cZ9oKnS2jRFerzcf6c4mJmtoAONC6e3tpbW1lbKyMmpqaqwWzNmCSbq6rtPb20tLSwtlZWWcdNJJuN3uad3WZKRrn1hRV1fHli1b6O3tzdnRLFX09nZ3pEmSBN/5Ds4bbsAFOHbsYD6Jucuenp5x+WJTApdMXtncCCbCXDe8mUgulmkrtH0qdCoyzqUbLU+6s4BMGhrs02wrKiqs7rGurq4pWSbmEk2IokhPTw/9/f3MmzePDRs2TDvZ2reVquBnkm1PTw81NTUJEytyIUnzeOu6Pi53OhXSjX3hC1ZUKfh8OG+9lejll2e1Dvmhh1K+ny53aY6K7+3tJRQK8dJLLwFYN+eFY2MUxCVXU8FMk+5MWzvmotFNboU2kW4qtBmctLa2pszTJ2N0dDRPujMJk2xHRkbw+/0sWrQoZfdYR0cH3d3dLFiwYFz3mCzL1p02W8iynNWFrWka3d3dDAwMUFFRkXKM+WTI9ouaTLrJZLtly5ZxX5ypFAmdN96IeOiQ8fqGG4ju2DErExjSQX7oIYTubgC0hoZJ12V3rDLnxS1evDiBFLo/+lGKHniAaDRK74knUnjwYEKElulAzLmeXpjOxoh0rdCDg4P09/fj8XgS8vTpWqF9Ph+LFi2aln2aScw50k1uaFAUBZ/PR3V1tfU3pi9CX1/fhAMWp6K3Tc7NpoOmaXR2dlpjzBctWmRpg3PZXjYRlqlEMG8+XV1d1lieiR4NjyWvhNgXvoDz1lsRGhsBEqLd6YgUXVdcAUBk166E9+03uHGkcM01ACxS1YQJB/aBmMn64uTztvC3v8W9cCH6v/3blD9DKsxWTncmoWkaBQUFLFiwIOF9+3h4r9dLZ2cnl19+OSMjI6xYsYLe3l5OOOEEzjrrrJyvkY6ODi644AL6+voQBIGLL76Ybdu2TcfHmnukq6qq5WspiqLRThpPEYTDYdra2hgcHMxo9thUSHeyZVVVtYjOHmUfPXo0Z3vHXCLQoaEhuru7x6URpns7ANFvfhPXjTcar6chvTBVKOefb6UXlPPPn/b1S5KUdiCmmbe0a11dLheLH32Usv37cXo8SOEwwl13Ebvoomnft9mIdGfD7CZVkJGqFfovf/kLl1xyCWeccQaSJPHMM89w9tln57xtWZa5+eab2bBhA2NjY2zcuJGzzjqL1atX57xOa91TXsMsw3RmMiHLMpFIhDfeeIPR0VGWLFky4WQEO6Yj0k2GPaqsrq4eF2XPlNOYHaqq0t7eTnt7OwUFBRNGtsk4Fq0do5dfjvPWW63XychWg2vCdcUVyI89BoD4D/8AQOivfwWmlnNN1Xig6zrCXXehDg8jDA/jCQYJNjczIMuMvP76uCLSse6neyya3QQCAU499VRWrVo15W0vWrTISlUUFxezatUqurq63p2ka4ff7+fo0aP4fD6WLFmSlVk3TG+ka8+XVldXs2XLlrR36ZmKdM3ourOzk5qaGo4//ni8Xm9WX46pmpgnE990RbrpCmhFu3YhtrQAb+WR7cgowh0dRYiniZzXX0/0qqtmxE/XIcs4PR7YsoXA/v2Ulpfj+N73KI/ri+0FvOS8ZVFRUcb5YnhnkK6iKFkVmUdHR6fVS9dEa2srr7zyyrQYmMMcJV2fz0dzczOKorBs2TICgUBCAj5TTEekm2qM+UQXoyRJOSkmJiLrZLI1CX94eDinjrRsSFLTNDo6Oujo6LByn6YlYTajzOX77wdA+fzns9pfAPHIEWNfVq7Majm9rg7d40EIBqGgIOvls4VdjfHmRz7CySefnFDRtyOdhaMsy1Z+2TzGqeoK7xTSfbslY36/n/PPP59bb7112lqM5xzpRiIRmpubWbZsmXWAj8S/dNliqu28nZ2dNDY2ZtUqK0mSNRkgG6TaV3veeNGiReOi61ytHTNZxi7DW7hwIRs2bACwiGJwcJDW1lYCgQAOh4NAIJBAFAkpovvvR4y38sr3358V8eoFBRD/zHoOemtt0yYAhK4utJUriV51lbGuJBewZL+IXGHpevftmzBqTWfhmM6oxul0JkTFsyFJm6psbjJkS7qxWGxa587FYjHOP/98Pv/5z/NP//RP07beOUe6brc7a+OXdMjlojTHmPf09DBv3rysfQmmI6erqiqdnZ10dnayaNGitOqMmTAx13Wd7u5uWltbqaqqsgqEZgdgcmGpo6MDTdMoLi7G7/fT0dGR0J1kaV+z/ILZoeZ4PUS//nWct9yC9NxzKGefTeyyyxI+pwnntdcivfmm9RqmTr65Il2+OBqN4vf7CQQC1jHeu3fvuKYDc5jqVGEWBmcS2UTT012s1XWdL33pS6xatYqvf/3r07ruOUe66S6Ymb6z2311ly5dSlFREdFoNOtHuKnkdM123cnI1r7MdJmY2zvo5s+fzymnnJKR7M0sfM6bNy9hsoXdK6HrrLMo/NWviMVi9B9/PEWHDyekKSb6jIGvfAXXz38OQOyyy3DccYfx+uKLM/q8utOJVluLePQojttuSyDe2bAtnA7YvRHMdtxAIMDJJ5+ctgMsOUWRrZn8sZhegOk7Z3//+9+57777OPHEE1m/fj0AO3fu5Nxzz53yuucc6cJ4YjBzs7mO6ZiIsMPhMM3Nzfh8vgRf3b6+vpwaK3IhXXPsu9muOxnZmpiOSNf0qWhqaqKsrCyjWXB2pEtX2LWvCxYsgHgBbFHcW9fv99PX10dTU5MVVZkkXFRUlDD80yRbxx13WEU1xx13ZEy8qWC/JlLZUaaCKx4RRW65JeftTicm6gAzUxSmztUc+WNvOCgqKkr7nTrWSDccDk9r5P3e9753xqSOc5J0k2E+3uZCuumaDuwm5suWLWPVqlUJxJxrES4b0rU3VrhcLurr66mrq8t4W1MhXV3XGRoa4ujRoxQVFbF+/fq0U44nQrbqhVTeuvb2XLuBSjgcJhqNUl5ezsJAAJfLhQNIdftMJyuLXXYZjttus16nQyZkK732mvXzsUK8qSCKYtqRP+YNb2BggJaWFitfnGwONBukm00xcK44jME7hHSnQ/plkq45xjwQCKQ1MYeZ1dvayXbhwoVs3ryZ7ng7azbINb0QjUbZt28fbrd7nOuZdPXVAKg7d2a8vqlGDPb2XLuBysGDBykvL0fXddo/9jH8fj9Vv/oVoigy+v73U9TdTVFREfN+/GPE+HilVLKyVDAjXfm++wBQ/vVfJ/x78bXXEHp60I+BNtRcj3e6G14sFrPyxV1dXQQCAeumV1paahHyRI5hueCd6KULc5R0k7/I9q60bGGSrt/vp7m5mVAoxPLly6moqJjwpM9EpGsn22SviHTmNRMhW9L1er0cOXKEUCjEpk2bxrmvSVdfjfj668YPV1+NunMn0re/DUDsW99Ku96ZekwzH58TpDzr1xtV7HhU3NXVxby770aPRBg9+WSUkREC/f0WSQiCkDLC1XWdkl//GnFkBAD5vvvSEq+2ZAlEowheLzpvf3phur10nU7nuJz8q6++yrJlyyy3toGBASvdluyLMJFJzUTI5roZHR3NR7qzialEupqmcfDgQXRdt8aYZ3KBTGeka5dfpTLmMZeLRqNZbStT0h0dHaWxsRFBEGhoaODIkSMZ2V1K3/42gq2qr+v6uMfwdLrf6cp/Fjz4IA6Ph9gFF1jvORwOysvLKS8vp2DNGoRQCHSdssOH8W3caMmtzCaEwsJCFv3977jdbsRPfzrnabnqli0ZGevMNGZjYKSmaRQWFlJaWjouX2zKBn0+X4KZfHLxbqLjnG1h3Ov15kl3JpF8MnKJdEdHR2lqamJ0dJTFixezbNmyrJbPlejtRJgJ2aZaLpttTRQtjI2NWV4QDQ0NlJaWomma8aWNP1ZrtuhO3bkTbOkFM8rNBjOR/5T+/GeE7u5xaQPnf/0XQiiEoCjo0SiiolBUWJhwrtW4aQ2iyOjoKJ2vv25NGxl93/uo/+tfcbndiJ//POloLLZtG44f/tB6PRkS5Gg33wxA9Mors/zU6TFbpJtqG6IoWsVOOzIdPlpQUGA5+GWTM86nF2YZsixnTLrJY8yHhoZyKhDlGumaqZGOjg7a29sTtK4zsb1UCAQCHD16lEgkQkNDQ0LrpCAILHj0UYR4FJJMvvZcrvrd776VXvj2t1FVdVwRayYNbyoefhjPk08ier2Ivb3A+EKZeuaZSE8+CbEY6ubNROOuYiYkSWLeU08hxBtWqtvbiX384zQ2NuJwOBj6+MeNfGbJMbwGAAAgAElEQVSSr675L5XUyrV9OwCRm25Kud+6rrPkwQdx/u1vVjed8+abp414Z4N0sx0hP9Hw0eR8saZpOBwOIpEIfX19GY37mStTI2COkm6qSDcYDE64zMjICE1NTYiiyIoVK6yT7/P5co5YsyUT01M3EAgQDoczIlv79qbiiQCGIqOpqYlgMMiKFStSplLG/bxvHwCSzweAmpT/VL/7XeNFmmOYTLpmVJtJesH5/e8DEM1gppnQ34/44osJ70X/4z9wf+pT6GVlRs51bAzHT35C7JJL0ioW5N/9DvHwYZYcPWp4Hfz4x9bv7Nri5Efnog9/mKKiImq+9jWEw4cN7ez27QnE67j9dgDCF1+cUmExXZgN0p0upMoXm5Ow29raCIVClloF0g8fHR0dZcWKFW/Xx8gKc5J0k5HuUd88ec3NzTgcDlauXDlOJpNNlJwrTLJta2ujqqqKwsJCGrLM/U0l0rVrjVesWDFpkbD3Qx+iYfdu44flywEQXn0VPS4SH7dv27cjahqhFIqGZNL1nHkmAKEnn7TeS545ZpKtePgwQlsb7gMHCJv7E4fzxhtxDA4SOv98HPfei7ZwIXp9Pa4vfQmAyN13A6CecQa4XIjNzQijo0jPPYf07LPoca8OsyFC/t3vkH/3OxBF5N27KQ4EjBvrlVeivv/9ACgf+1iitjgOsxAbCAQIhkKogYAxmmZoiP6mJoqKiqjavRuho8MoTP30p2geD9Err5yz6YWZhDmp25zXZyLd8NFrr70WRVHo6upi4cKFnHjiiTl5sZj44he/yB//+Eeqqqp43SwcTyPeEaSbnNM1NaZNTU243W5WrVqVtjAky3JOXgiZwE62lZWVVhfXwMBA1uvKJdKNRCKEw2Fefvllli9fPk5rnCmEZ5+FSAShqQnpttsSol1p+3bE115D13VcO3ZAPFUT/c//NJa1ka7nzDMR4/PFPGeeSejJJ42ZY3HfBZN8re22tRlmND5fgieD88YbEQ8fxjM2hlRcTGzrVuQ//Qnx0CF0SUJsb8f1hS8QueceYpdcAqKIMDgIgF5SgvjmmxbpmlD+8R8RDx+2pGXZQJblt1pz77wT1/bt6EDsu9+lNK57HRoaQo6fdyU+/2toaIiiyy4zIuqst5oec510IXVjRLrhow8++CBXXnklFRUVPPLII/z4xz/moYceyrk7bevWrVx22WVcYCvOTifmJOkmH0wz0s1ljPlMRLqaptHT00NbWxsVFRXjWmaTBx5K8UdQNZ4LTIVsIt1oNEprayuDg4NIksSpp56a9QWo/vd/AyD+6lcIZWWQ9ISQCvIf/wglJeh1dTj/8z8t4jUhdnUhBALoKc6J/MtfgttN7Gtfs9IJ7gMHwOdDfd/7Jtxu7CtfQdy719jG0BAUFiLt24fnE59A+eAHiV16KeIbbyC0tqKtXo22ejXSgQOoJ56YkF6IfvObOG+8EW3lSkYbG3G53Tje/35r5I/8hz+gfOxjkx6HyE034bz+ekpvvRXPVVcZ2uLvfhfH7bcb05//+Z/pampi3sgIHR0dRKPRhOnE5qNzrs0HM924MBum9Nl0o5WUlKCqKp/73OdYs2bNlLd9+umn0xoPDmYCc5J0kyHLMoFAgD179lBcXJzVGPOpyM2SyTOZbNPNQZMkKetoJJNIV1EUWltb6evrY8mSJWzZsoU9e/ZM3RD7Qx9C6OsDUuR0b7oJ9YoriLz0EkJpKS6fD1pbEeLyNkEQWPjjH+OcNw9140akZ54BhwNt0yZcV11F5PrrcX/ucwiahu5yIQwPW+sO796d0vLRJMfw4CChbdso5a10QsGmTRAMgiQhvvACzj17jDTF2BhCKIR48CDaqlUoH/4wyic+Me6zinv2ANDxgx9QUVFB5bPPjvsb+Te/MfbJ5jzluPNOAGJf/nLa4xj793839jEUwu12J+QgTb2rqS02C0r2wl2mBufvFFvHbNULM+GlOxOYk6RrXnR2E5ZoNMrJJ5+ctRJhurrZTLKdP3/+pEMnJUmylpNuugkh/kg7UcQ7UaRrN1BPtplMNUY8G2if+lTa30WjUVpaWhj61KdYU1ODq6sL4be/RQkEaF+/ns69e5EkiSWRCNqTTyJ5vaDroGmIb7wBgOuqqxC8XhAEBK8Xx+23IzQ3E4kXsMSjR1Nv+5vfpOfwYRZCglxLW7YMsaUFobMTIRIBXUd64gli8WMqtLQgvfCC4Z1rI13Hj36E9Je/IAwPo61fT9V996FdcQXKxz6G/Ic/GMc5gygXDCN0+xRi0y7ShFn5d8aLiNGvfz3tdGIzh2lOybVri+0qCntB9p1AuqqqZqWXzpPuDMNuL1heXs6GDRt4+eWXc5J+mZrAXCBJEt3d3XR3d2dEtvblss3PplrGbl6ezkDdjJCn80toup319vayZMkSjjvuOMSnnzZysDU1yPX1LC0sZNEppzAwMEDnxRdTeO+9uPfuRXc6GV2zhpK2NmRJQnzzTSguRigsNKJpXUf+y1+QTjsNMIzGwWi+SOV/UPTTnyJ2dAAG+YZ/+Uvj/QULDEWFw4EgSYivvory4Q8jdXeDriMeOoTzuuuIxodMmhD6+xHfeAP9tNOsm7udbOXf/OatdEM84hWGhiyjHc/HPw6hEGJXFwDSX/4CcdI1HdAKhoepGRlBHB01PtsttxBNYR9oN6yxF4ZMbXEqj4SioiIrSpwp8p2tSDfTp1Uw6he5fP/fDsxJ0jVdkrJ1vEqFXCJdXdfp6elhZGQEh8OR9Tj1BG/c7dszyuna3brs7cKpzMvtWHnjjTjnzUO7667J9yuuGuCMM1L+3j4lora2NtFLOMUjryAIVO3aRYHfz7y+PqitRV29Gv1b30Lcvh39yBEG166l8NlncQWDyC4XoqYhhMMIw8MIkQhCIIB2/PHWOh3/+78ACR1odrji7cja4sWInZ0giigf+Qja8uWo556L9PLL45Zx/uAHSC+8gDAwYOh1u7sp2r8ffv5zyGICrPTMM4jNzcZnNN+M5z/tDmhyYyPCFORNqQZi2j11u7u7GRsb46W4ttgus0qnLc4Gx1p64ViaXp0J5iTpyrKcUnKVi6duNp4GJtm2trYyb948Kisrqaury2mcup3oJyJbE6YKoLOzk7a2NquDrSDu7xl74olxy8gXX0zR4cOILhfixRejxCOtlPv0/e8jHDoEwNKODnRb8c3+ZLFw4cKUJK984xvEYjHcpgTK/EwpmiMkScLlciG63Szavx8hFEJfsADd7wefj2h5OV3nn4/zyBHKDhwgsGABvi99iQU//SlyX58xBfp//xc2b0atqUFobDSKYtu2WaSr19Wh1tVBLIbQ0QG1tcb7NTWIDz+MXl09Lsq1jrXfj+zzQVvbOItI5Z/+6a0IN65GiH3lK3g++lEIBNDLy0FRLOMb9fTTx60/cvrp9J59NlV/+pNxrKbBJNvuqWumIKqrq8dpi7u6uhJsHO354kwLV8fiUEqYPi/dz372szz99NMMDg5SW1vLtddey5fiUsTpwJwkXZg+T91MTpQ9dzxv3jwrwj58+HDOpjfZpBfM7ZttlKYawnHWWYjPPAOA46yzUhIvgE5qu0MT4u23G00QNoWCeWz7+/tpamrKyrjcjrGrrqKjo4P5/+//AW9JySLXXw+A+/OfNwZDhsMIHg/6woXIVVUsGRtDHBpCX7QI0eViNC6AH+3vR9c0FEnC3ddHNBQivHo1Um0tzv/6L/TSUqL/8R8W+Ua+9z0jX7t3L9xzD0Jrq0Hw9s4oWUZ9z3tAUZBeeAHd7Ubv6sL1xhvotijbOl6trYgvvIAe15A6/ud/0NasQYwfa23NGkuhEY3vR+zii630gv+f/xkh3kE3E9A0zSKscb7FcaQbE+92uxOi4lTOYceal26utq7p8OCDD07bulJhzpJuMqbiqZsO6cjWxEyPU9d1nf7+fpqbmykrK6OgoICVWQxPVO64g+AnP4lUUoIwQZQLcYVCXDnQ8cEPIg8O0tzcTHFxMac98ACyJKHEC1YTIZoUtVtG4HGydfz0pwDE/u3fAAjffz/O730Pad8+tPp6Y1/q6622Xqmx0ciB1tQgNzQgdXejnnQSY5/+NEO3304sFmNwaAjP00+Dw4Esy2jf+haqy4W+YweeH/0I8ehRpJdeQnruObRVq4zCnaqOmxQR/cY3rOYKNM34Z7upOm+5Belvf0MYGgKHA314GC0+Gy5yww2GTjn+OhXMiFkbGUEUxZwiXMePfmSs66tfTfs3mZBiurE/4XA4IV9sOoct/b//w+VyoWzdSiQSmfGpGtmQrs/nm7ahkbOBOUu66SLdXGFPTdjJtry8PG3ueKaMzHVdZ3BwkKamJoqLiznppJNwu90899xzCX8Xe+IJHGedZb0GkHbtAkCNewy0fec7LFmyhHSXpHj77QhxLwq9vp7hz34W/0sv0dXVxYknnmg8dk4hqrGfJ/cFFxhNCpqG9MwzhONyMDMadH/xi8ZnuewyXN/4hkFumoZeVITn/PNRzjgDbcUKKCmh+NFHaf/IR6jat4/ip59G6ulBKywkumwZStxM5fX9+6lpb2fh3/+OrutIgN7bC5WVEIkYjRnRKGJjI4ClJrBD+r//QxgdJbpjh0G4g4MIo6PopaXolZVoy5cT+8pXsjomMz1aKtcCmiAIeDwePB5Pgm+xqqrw8suEIxH6h4cZGhoiFovh8/kSzGqmoi1ORjY53blkdgNzmHSTMRVPXXult6+vj5aWFsrKytiwYQNut3vC5aY70jWnNRQUFGSkN06XUjCRaSdbLBajt6ODnqYmCgsLWb16NS6XC3nbNoT9+wGQ40Ul5Yc/RIwX5rSLLppwvSbpOu64w6jyxwle03UccVlY7NJLkR98EG3VKmM78cc79dRTkZ54wug0mz8fKe4DodTUWOuWurqMScKShDg4iEsQcNbXo551FpsGBlB+8AOc552HdOAASnEx4eJiw293bAxJFNGfeAKhuhpnRwdiIIB2wgkAREdHcUajUFJiKB1uuAH1Pe9BPHIE6a9/BVG0lBLAhFGuM+5PEY2bA2malhPpOuJRu/k6XbQ73aoF9y9+gdDXhwco+/vfGRkZwfepT1FZWZnSrCYXbXEysjHUmUtTI2AOk266rrRcYEq/Ojs7KSsrsyLLySDLcs5z0pJvECMjIxyNm6yk66SbTHMr7dqFEI/apF27UK+4YlLS9V94IaM7dxKTJIq2b2djeTkvv/xy+mXefBPxrrsQ4h078mc+g755M+rXvpbyz+2RrrZ2LUJXF+g6gs+H/PDDKB//eMrlIj/4AQAF69eD04l6+umIr72GdsIJbzU1xB26TOi1tUb6YGzMes91441I4TBiSQlOQI7npNWTTkJRFIJ1dUSjUZz9/QwcdxxFfX04nU58GzYghUIUtbcjNjdDIIB26qmIL7yAMDaGvmCBJfVy7diB+PrriC0tiHv2EHr6aWv7zu9+FzFeoDTJV//3f6fqF7/AMX/+hCOCcsVM63Q1TUOSpLRmNem0xfaoeKL5a5BdUWwuOYzBHCbdZOQS6eq6Tl9fHz6fD6fTmTHZmphKpBuJRADjLn306FFEUeT4448fZ8hjRy6a23SkG4lEaGpqMkxwtm1j0dNPIzz5JNonP5ngoKb88IcpI1wwHMiEWAzhoYdwPv00wnvfi5JCYqXrulVI0i65BPnXv0bs7obiYsTmZlxXXknk5puRH3wQ8ZVX4PjjLR8Gva4ORkaQnn4aIRYDUcRx773Gik87jbGvfAW32430wguoZ5xB7MtfRnr0UeRHH8Xx5JNgFv4EAW3ZMjSzTbS8HBnwXHklxbt3Iw0OUjJ/Phw4gKKqvLF9O/1OJ6u+8Q1Km5tRh4YQDx9GDwbRystBVZGeegpnMIi4Zw/iwAAEAoiSZAyxNFNfKc6Vrut4Dh1CnKRFPRmxr341o5zudJNu7MILcZgTly+8kKGWForSPPrnqi02/01m4ZgKc6kxAt5BpJuNh4K9QFVaWkpFRQVLlizJinDNbeaa0w0Gg7zyyitompZgNTnZchNJadQrrhiX000m3VgsRktLC4ODg9bATemhhyzBv/jrXyM0NCQsYy+gmekE8a67EI4ehc5O8PkQdR35mWfGkW7yF8hx332GhKuoCKGzEykWQ9240fhlURF6fT1CdzfioUNoa9eiLVyIODqKODaGXlKCMDSE/MgjKB/5yFufads2YoDjJz/B8ZOfgCwjvfwy4sgIuqaBx4O2fDnamjVErrvOWk7evRt5927Ezk4Erxfx5ZfR16/HgTGifPXevRQGg4gOB1IggCYIaEBUlonKMpokIb/yCvrmzRQ9/zyiJKHHUyQmot/+9rj0QuXnP4/c3o7odOK66CLUs89G+fSnU57TZExEtiZmItKNXXih9ToX9cJk2mLTFCgYDFoFvdbW1oy0xaOjo/lC2mwgF0/dZLI1I9upSL+yjXQDgQAtLS0Eg0HWrVuX1R06E6mZmmTSbZKuqqq0tbXR09PD4sWL2bJlS8ovpvD44yzYvx/9G9+YdH+0D34Q8Y9/RAC01atR4uqDZJhRs/Tcc4ZCorAQXZahogLl7LOJxqdRJHyOD3wAraoK4fXXweVCq621FA0MDhqdbPGutXTQS0shGEQvKCAUz30nj1IXX38d6e9/B1E0HMYOHkRbsYLa/n6jrTkYNBomJAmhrg4xGsVdX4/L7YYjR4zxNMBYQQFSKIS3vBzn4CDznnkGeWQE8bnnCD/+OGBIywAikx7ZqWGutAHbtcXz58+33o9EIhw4cAC3252gLTbtHpO1xT6fj8WLF095f2YLc5Z0kzFR1Gkn25KSknHjxHONWLNZzm4gvmjRIsbGxrJ+JMplOKUgCAwMDNDU1ER1dXXKVmHtk59EMhsASktxd3biuvNOmEDSZEa82kUXId18sxG1XHbZuJE2qSZHCD4fxFMY4uHDOHfuJHr11eN8DsQXXkD94AdRP/hBpGefRTQLVyMjOO6/n3X33YeyYQPC5s3oFRVG7hUMRcG//AvyH/6Atm4dkeuuw3XNNUjPPYcen+flvPZa1NNOQ3rlFcT2drTqaqORor8fVqygdO9eRLcbJAndvMHPn4/Q3o7u8yH09hrtxZJEUWEhyoUXopx3HiWKgufcc5H7+owmiUOHGPjqVw2fjcpKXC4XPd//PnVXXYVcWEgkqVPQFW8ZNnXMuWCmdbQzvX5N03C5XCxcuDDh/VTa4quvvppgMMjKlSvxeDysXbs2K1llMh599FG2bduGqqpcdNFF7MhgcnS2mLOkm8mcNNPqsampKSXZmkjuEMsUmUS64XCYpqYmxsbGqK+vp6KigkAggNfrzWl7mTZVmN1zbW1tlJSUsHnz5gl1j/pJJwEgNDdPOmJHfOABALTPfQ4A9corjf1KMTjTUi/ce69lBi6++ablX2uHGYXqlZXWhAfxhReM9+KNCILfb0nJJF03iPTQIbQ1a6yhkOIbb6Cecw7aunUAuK65xrB2PHoUoakJbcsWwIg8xY4OKCw0VAHxm6D4yiuMbdmCu7XVSLsUFqLHGwu05cshGkWvrTWidSAcb00G40YsfPjDCE1NhuFObS11PT2oqkpYVYlFozj+53/Yc9lluN1uig4ftiK3ip07EQ8eNPY57sCWC+ZKpJvt+lNpi//85z/z1a9+lVWrVnHkyBEef/xx7ow7vuWy3UsvvZQnnniC2tpaTjnlFM477zxWr16d82dJhTlLusmwR512X93i4uK0ZJtq2Vy3mYxIJEJzczNer5fly5ezevVq60aRS8Sa6XL2G01ZWRlLly5FkqRJheZa3HZQvP12ogUFqF/4AqlcLcQHHjBysoC8dStg2Cq2tbXh8/koLi62SMTpdI4j8LBZBAOc8UkT0auvxnnttUhvvmmoIlwu1E2bcG/dCkNDRK+5xrJMlO+6y3ARixci0TTDylFVDUIE9PijZuS663DefDPS449DQQF6RQVCOIy6ejXR73wH1xVXoMf9EARBQBgZQV+wACEUovzZZ42o3e02GiE0DeXcc42Ie3QUhobQq6sJxwt+7vixCN97r9X0ID39NOpppxmtwUBBczN4PCilpSxfvtySXPn9fnp6etAHBij0+Sjs7ETr7mbg61+nqKiI4gcfNEbFxydrTIa5TrrZNEaIokgkEuEjH/kIa9eundJ29+7dy4oVK1gev44+85nP8PDDD+dJNx0cDgfRaNRKIxQVFbFu3bqMnIdkWbbUBNkgVaRr2R0ODbFs2TKOP/74lPK2XFUPE0W6w8PDNDY2UlhYaN1ourq6slZ12M11kiE+8gj4fOjve58lD9q3Zw+LFi2irq6OYDCI87rr8EejHL3gApxOJ/N/8xuC0ShOhwP5nntQ4uQRvfpqnLfcQsEHPgCxGHp19VvbOXgQJAkEAecdd6AeOEDsy19Gj+uWXTt2oAeDCLpuENlZZyEEAsi/+IXR5rt4MY5du9BXrDDysgcPotfUoHs8SPv24fjhD4ns2oXriiuQH3kESkthcBChowO9uBhRURCjUYjLmvT6esTGRoTubrT6eoTOTmtf3Vu3Wvpj99at6G432oYNhH7/e+tvHPfeizAwgLp5M8PnnIMzPpLGHrm5SkuRYjEQBGIFBVR9+tP0fuADBNvbjdx8by/hCy6wbmrpiG+mmy9Ska7jJz8BMCZ1TBHZkC5Mn3qhq6uLurirHUBtbS0vxJ+0phNzlnTtF5U5C21sbIz+/v6sTMwh90jXvg+xWIzW1lb6+/tZunQpDQ0N6fW0OUrN0kW6Pp+PxsZGZFnmhBNOSBhNZJenZYRFi1ABIQXpSjfdBB4PQk8P6i9/ycCKFbjcbt6zdy/atm1Eo1Eq/vu/Efv7AVj4+ON4v/51Is3NSOEw3upqxlpb6XvpJQoLC1myezfuhx823MBkGQ1Qzj0XvbISae9etGXLkF59FWFwELGlBcdddxGL55KdP/oRelubQZCCgOPBB9HWrDEI9vBhtOJihEAA8dFHDfKOxRAGBtBra41mh0OH3hqXHggYxjwul5EicbvxbtjAvHAYqb0dfdky9IoKpP37wetFePFFtBUrUE8+2ZirZoPQ1IRQVQUtLePMctTTTye2dStaayuiKOKOR+9h2+Owunkz0muv4cDQFNfW1hoj5DUN4bHHCB08SOPVV+P3+8c1IhQVFVkKnGOx4y1TZEu6o6OjecnYbMHeLmuac+QyrsPhcOTcWKHrOk1NTfT29rJ48eJEu8M0mEp6wR6B+v1+jh49iqIoNDQ0pJSdmctIcUct1SaZGve3v/sddHfjGBhA/tOfIP7IbEc0GiXq9yPHYlRUVCDLMrokoZH6i15++eV4vV4KFi2icO9eqt7zHhauW0cgEMC9f79BuNEoeiyG1tqKz+slcOmlFH/pSxQ+8ABiezt6UkHF9a1vIYyNobpcVuFO6O5GCocN0oxEEI4cMcby9Pcb73k86OXlRo42fn7cX/4yuiiinH8+8p//jO7xGKSu6wS7uykpKECIPynJjzxiFP/cbrRly1A/8AEcd90FLhfBPXus9IJeXGxN2TDhuPfet+bA3XsvtY89huxwIMWHHrq//GXCd95JJC73c8UVKPqSJQg+H/rSpTjuvx9CIUr27GHDpZcSfPLJcY0IPT09hMNhQqEQh2254omi4lxhP9eOn/zEivTNXPxUkG36IhgMTouXbk1NDR3x1BlAZ2cnNfHux+nEnCVdXdfZt2+fVbEsKCgY502QKXIppJkG4oFAAEmSUqoC0iHXKMSMkEOhEE1NTQQCARoaGhI6glJtK1vD9FSFNL/fT+MHP8iiri4qzjsP+ZprkP77v9EB9WtfQ961C0FViV5zDc44sevFxUi9vZS+8QZie7tRBNu/n8Lt2/GceCLCiSciPv+8sU2Hg9h55zG6YweF//Vf+CMR9n/2szjWr2fTBRcgiiLBa66h+LOfNfK5koTudCLEYqjvfz/ioUOG1aLbjb50qeGR0NUFsgw+HyiKEb07HIZ6QtOQHnsM5ZxzjH2trEQYGUH+05/Qq6spGhlBiE8/Fg8fRi8tRRgaQgiFYNUqHA88gDA0ZMyD27UL7cQTjVFAJ5xgjKv3+60oV3riCQSvF3XTJmNbTOz6lozYF76A9MILSE89ZdxA/H4r6k/ViLB3714WLFhg5Yoniopn2rgmFyiKkrFm3rxOpyPyPuWUU2hsbKSlpYWamhp2797NA/Gi8XRizpKuKIqcdNJJ41oJc8lnZZNesBuIV1dXG4/JS5bMysWraRpdXV20t7dTX19PZWXlpNstu+EGHIGAZSUoXXNN2mhX+8d/RPzd71AkidiHPkQphvri6NGjBAIBjjvuOMpvvdX6++TWX/G115AfesjyqXX88IdoS5Yg7tuHEAhAUspH+utfjc4tUUR3OBBOO426n//c6O4Cam67zej2GhpCEwQ8//mfxFQVKRRCB+S4Llt88UX05cvRq6oMyZeuGxrfoSFjQ4JgbCcYhOJi9FjMOG66jtjXB+EwekkJ+P1GJFxQwPCppyJfcQWFN9xgGKp7vTA6aoxzP3AAwes19LsdHUhPPUXoD39A/v3vrSYTE85du4x8cX8/8m9/i3rOObR973vMmzcP1ze/CSSmFwAiu3bhvPVWq6XbeeutaMcfbxjz+P2ocZOjdFiyezdVdXWU2TTb6aJiWZYtzWtxcXFOpjWxSy55W3O6MD3pFFmWue222/jQhz6Eqqp88Ytf5IS4F8d0Ys6SLoDT6ZwWT91MSNc+Tn3BggWWBKuvr8/qRZ8pKIpCS0sLXV1dlJeXs3bt2owvskz+TrznHqTdu9FXrUK59VaCbW0Qi3HkyBEGBwepr6/nhBNOSLsu6dZbER5/HH1sDOmmm5CeeorIj39sFLZ+9StUUUT56EcR33wTbe1ag1RuuslQC5SUGIQcfzwUX3oJ6eBBtAULEHw+Y7qvosAk4JEAACAASURBVCACLq8X9X3vMwZNdnejCQK6wwGDg0RVldiiRThEEdHpRKiuRhYEo/A2NobY1IReXIy+dKlB0oJgNDy0taE1NCCMjBi5aEVBUxScy5bh+vWvkV58EaGnxxi06XKhl5UZrcVFRQaJx2Komzcb5+m885B//3vQdUsPbJ0Hnw8BQ85WffPNRHbuHEe2kyH45JNGSgOs3HYynLt2UdDZiRgO49y1i2iceNO158ZisZSmNea0CVONYnaEpZMSTgfZmsjGwDzbAZaT4dxzz+Xc+GCAmcKcJt1k5OqpO1Fhyz4toqKigk2bNiWs3yTsmSBdVVVpb2+nu7uburo6jjvuOEKhkEV+UjxiVdNMQAAIfvvbdLa3Uxn3r1WSolzxnnuQ7rzTyEN2dCBt28bw1q14vV5WrFiRtnMtGUJXF4TDIIpIzz+P+5OfRDpwACESQQL4y1/QCwuJ7NqF+zOfMQiuogJ9cBC9sJDo9u3IDz+M4PUa+tdwGG3lSoSRESNtIMto69ej/Mu/IL3wAvI99xjR6/z5CA0NeA4dwqlpRBYtIqIoCK2tdK9bR7S6mnmHD1PW14coy4ZnbziM2NNjRK6qaqQjTN20qiL4/cx7/nkcIyNoy5cjjo0ZRbhlyyAUQt24EX3VKuSf/Qy9tpboNdfguvJK43gePIh45EjCQE2toQGhqcko1JGZg1b08stxxp8qopdfPu73pv9ELCnvnu3oGofDQXl5+biBmMFgkEAgkLIjLBaLMTo6Oq1WjnZk832aay3AMMdJd7o8dVNFcKYZTnNzM/PmzZtwnHqugy3TpULsaQRzBpokSfT391v5Wem666zxOtJ116UlXrOQpr3vfZPuj6qq9MXTEOteeomKtja0FMW0cctdfjnC88+jd3YauVHbGHVkGV1VobAQbcUKHD/7mUFwfj8UFaH8678CIASDRtGspAS9shK9rAzlYx9DLy9H+r//g/nzrRHm6DrK1q2EXnwRp9OJ8tBDhhH6U09R+Pzz6OXlaHV1uGtqCFx5pVFw9Pspuv12otEongULKPN6cUajRpR99CjayScbkjBZhuJiHKOjqOvXI7W0GOkPlwtt+XKjkLZpE8LgINrGjYgHDlC0cCG6KEJZmfHZBQH8fuQ//xnlnHOMaPOKK6xus44LL6Q6ft7tRjLJSEW2sYsuIrkwZyfe8LZthLq70WpqrCg3WwiCYE2bSI6KfT4fw8PD46Jie654qjPYsjUwn0sOYzDHSTcZU/HUNWFvLigtLZ3UU3eqRub2i8tunp4qqs5lirAoitTceqsl4+Kqq1BtnU7q1q2M+f0I99+PunIlpXfeSWxwMOvPE9u927BIvPtu4+eLLsK5cyeOu+9GURT0z32O2LZtOH72M9RNm5D27kUvKzMe1wH5l7+EWAyxsxN1zRrLq1bo7zekWqKI0NVF7CtfQdqzB4DAunXoBQVIGKQtHT4MwSCCqqJ94ANEd+zAAW9Fcjt34v7nf0bw+dDKyhB7e1GdTvrmz+fIP/4jpz30kJGSkGVUhwNR140cra6jvuc9aCtWGGOFRkaMY3vggDH1NxRCwCiQ6YWFRtohLlOT//QnKCgg+h//QeT663HcdhtVv/gF4o4dOH7+87fI8+c/T0m82SJXr95M4HA4rPzvqrixjxkVJ89gM3PFdp+ETKPXbEl3Lnnpwhwn3en01LXLz+zNBZNhqkbmsixb2z569CilpaVpJ1XYpWbqNddklF6w2zQmw+fzceTIEVz/8A80XHQRxfHPW/zrX6N3dRk61/hjbCYRL7yVa3TcdpsRtRYUoEUiCPEvXOyLX8Txs5+hLF1qvE6T13Tcdpvx/29+g9DTg15YiPzHPyJ2dCC+8QYAUnU1ozt2sOB730PcswdiMSPlEIkgHjqEe+tW1FNPfSt1EAwaSoSFCxE9HoSiIqSSEqrq66l86inG9u/Hc9FFSK+/Tqi4GO691zjeRUVw6BCiz4ewapUxamfHDmPfIxGja01RoLTU+Ew//alRmHM4wNbu7bjtNsSmJtzDw3j+53/ANp0hG8S2biVdekHTNLr/5V+ojLdATzeS5Vz2qDh5BpvZbZdtVJyNZCxPum8zco10h4eHCYVCdHV1pTUQT4epRrojIyM0NjbidrtZt27dhE0duUa6LZdcwoL4BFv1+usJBAI0NjaiqmpKD18BY7LDhPsfJ0V1EhNu9cwzGe7vZ77tvVh8LA9gtfcCSI8/jrZxo0Fo8fULcdtFwmHLnFzo6QHAMzr6lmNXaamhgIhEQNMQDxwwiPLll6GsDEIhY8KEokB3tzF80lbsEgC32427vR0hFsMVDFJQXIxWW0vM6UQbHibS04PQ2IgOOH71K1SMcyIoCtq6dcQuuIDY1q3oJSU4br7Z2G9RRP7d75B//3tDz7x+Pe62NhyaRuy++yZML0yEZLK1jvcxYnaTKlecbjJxclScjQIpn16YZaSKdLMhXbOTyxxHffzxx6eMMCdCrpGupmkcOHAAh8PBqlWrJjQvN5HcVDFRhGtfRtM01OuvJxKJcPSNNxgbG6OhoSHBTs+O0Gc/C/fcg15VlTLClb/0JYS+PvSGBqTbbktJvELc/EZbsYKg00nJZZch33ef8bt4h1zsootwxbvCIj/8YQIBxy67zJhZputGpK7r6DU1KB/9KM4f/ADcbkInn0zRvfcS/da3cPzkJzjjml9jhREjFRAIoMdi6LYuPZxOtJUrUc88E4Bo3ElKevZZ1C1bkJ96Ck9PD6xahTg0hGPJErSzz6bg7rsRolEjf+v1gqahlJQghEKEvF76W1pw7dxJ+ItfZOnddyN5vUYkPDIC8+ejrV6NVl9PtLUV4cwzccZvLMnTI+SHHgJAOf/81Cd1AhzLvgsTTSa2R8XBYJC9e/dmFBXnSfdtRiaeumBUPI8ePYqu6xx33HGUlJTw6quvoihK1qSbbaQbDAZpbGxkdHSUFStWJPR6T4ZcI11VVWlsbGRgYGCc+U4qCILA0Mc/TmWSRZ5lrRgvhAmNjehpPHRNxC67jI59+6i97z7EtjZjMoSmoZ9wgiV/mgjqmWcivvoqhMPEPvc5w5B87VoErxdPUxOi34+OIVmSH3nEcB2TZXA60cNh9IoKY4RPMAjd3UYKQlEQgkGLbM1mDvX97ze61cxjIwgIkQji0BDCiy8iKIqxvCiir1+P+OqryKKIXl1NgcfDwoULER97jLHRUbyLF1MQCIAkISmKEWm/9BKOv/2NQk1D/u1vERUF7bjjrKnEJtmaWl/5oYeyJt5jmXTTITkq9vv9bNy4cdKoeGhoiMHBQZYsWTKt+zPTeEeR7mQEaG+bXbFiRcIdcqYm+5qwWzw2NDTgcrmy8oeA7NuHNU2jvb2dhT/7GZWVldR///sZfSEnm6um///sfXl4VOXZ/n3OnJnJNtn3hewTElAgIQit+KF16Q+1Klb93OinbdWvH4pFcKNaaEUUFMViK+BWrMul1kptFUSrIJUECigoZCP7nsky+3aW3x9n3sOZyUwyk8xMIOS+Li5CmDnznpk593ne57mf+/nBD4Cvvwbi471GuapNm9wn7PpQTjDvvusmrbJ7jHh3rFgB1aZNoFpawM2bB+fdd0O5Ywe4OXPAfPwxaIUC9h/+UJx9tnOnpHygXERMZrEpjh4VO9MAsZNNFvWqnnhC1P3qdGK+Nz4e/IwZ6J41C9nHjokdbomJYmTvIgVBoRBtJGkaVHs7hPh4UCyLyBdeAGW1IuLQIQhFRaI7WXQ0eIUCvM0G1umEJTYWDoaB0mQCA0Cw2UCzbMBSL184G0lXDlII9Ccq3rp1K/bv3w8A2L17N2bPno0HH3xw3IXE9957D2vWrMHJkydx8OBBzJ07d1zH88RZTbr+eOoCYnTZ0NAAm82GoqIir22z47F3HMlQRu46VlhYKEWZAwMDAacl/CV4oi1uamrCzL//HVxPDxJ4HsLGjeBcXVAjwZcInr/pJgBixEvFx0sevE6nE6dOnUJvb6/of/HVV1APDQHFxaILGAD29tvBvPEG+NxcUHY7mPffFzW4rp2JwksLt3LzZggKBayffCL9zrl0KaLOOw/U4CDYWbMgeMj4+PPOA9XUJHajkUKozQZeqwV3ySVgPv0UXEUF7C+84PXcuR/9CHRVFWxKJdiSEtCu4Zec7MYhREWB+fRTCNHRoiWk0SgaAel0ItErFEBDg+gZoVKBTk4GUlNBA+BefBED69YhLS0Ndp0OdrsdTRdfDPuhQ2Dy86HRaJA+NIQItRqK664bZgo/Gs520h1NuSCPirds2YIHH3wQ119/PbKyslBfXx8U5cbMmTPxwQcf4O677x73sbzhrCZdT3gSJ/EoMJlMKCoqQlJSks8PJdiRLsuyaGlpQXd3N/Ly8qDVat1eeyy5YH/Gqet0OtTX1yM+Ph6VlZWI3L8fgQrARnsd/qabxD88j4FXXkF/fz+ib74Z8+bNg8PhgPPyy4E9ezAUF4dTCxfCarGgrq4OeSdOQKVSgf/tbyUFAuWyP+QuvNCvtUX87GegDAZQANSdnWAhRrlUZye4iy4SdcAWC4TEREClEr1/WRaUIEBITITl3/92O57jN79B5JIlQGQknK6LzPHUU+g4dAiZlZVQbdwo/m7VKrfnUf394hw3rVaUkDEMqMREsctOoxGj7thYgKJAnzoFbvZsOFavBgB0//jHSJ89GzF/+hNiYmIQ61Ia0M89B7vDAXN8PPpNJrQdPoy0zz+HWq2G46c/9WuK7mQnXU/o9XokJSWhtLRUkrGNF8E6ji+c1aTrK9K12WxobGyEXq8ftYWVIFhG5jzPo62tDe3t7cjOzvbpOjZW0h1V/qVWu6kguIcegqW1FcK0aX5FuaO9DnC6ccTwxhtIdjhQkpgIxdNPw1laClVkJOiBAaCiAurcXKTV1eFrrRZFb7wB1NfDwbIw3Hsvmu+4A9HR0Sg2m8XPbeNGty+jcvNmKUWh3Lz5tA0jIOaEnU7wEREwX3cdor/6CgDAfPGF2HThWjvV0yNOgzCZQPX0QLV5MxTV1bD95S/icbdtEw/nyl3TJ064WTECkHTEnvAcOKnauBFCWhroEycguGa/KXfsAOXya2A+/lh87OrVEjGSAhqZ8Es3N0MJICoqCs7ly5H97ruglEo4bDaYPv4YXRdfLE3RjYiIkEhYo9FI5jXhIN1A6x6BHj8QUj/bJgEDZznpeoKMeT5y5AgKCgpQWlrq93ZjvJGuIAjo7OxEc3Mz0tPTRx2PM55ONjksriiSZVmUlJR4bYlsufVWZI0yxFEOX85kiqeegtVqxdHFi0W1R2EhlDodFK++KjY2RESI7l4DA+BdqQcASP/iC3FqgutGEJOZiZz77wfPcWh8/30xR/ftt+A4TqpYZ5nNiOI40B4XoO3Pf0bUD34Aqr0dVlcRj73mGjHa7ekRi11y0uF5sT3ZYgHFMKAbGxFx223g5s8/PZ4doieC4sABQKXyKcfyBJnuC4iRsGrjRvAlJaejYpqG8t13pWYKz/fY7b09cgRUT480ckj+OLVaDWVCAjRaLQBI03JJbrO7uxs2mw0KhQI0TUOpVIasTfdM89IdK+leeuml6CaDTmVYt24drrnmmoCPFwgmBenKDcRpmsaCBQvG5DQ21ukRJpMJBw4cQFJSEiorK722C3t7nsPLTDF/YbfbcerUKRgMhhHlX2OBt/QC97vfwXDkCCAImKtQQPH448CMGaB/9SvRblCpBH3sGLhrrwV98CD43Fywt90mPvnpp8WZaHV1YHbvhmLfPlAcBwWA4ptugsVVDCHdTUajES033oi4bdvE1uRFixBz6hQ0Gg1Sd+wQc6gKBRRmM4DTU3aFmBgxj0tREDIzIeTmgj52TCTX6GgIHvl+EtVG3HknKIMBQmIimL/9DYq9e5GVnY2IP/4RQkaG+Brbt8P5y19C6Sr2OT1GzQOAYt8+8QcX6RITGGbPHnBz50rpBU84770XwptvgvnsM/BFRdKYdfbGG8G8+670MwFFUYiMjERkZCRSZFpjch3Y7XapIUEQBDfplUaj8ev76QtnWnrBZDK5mfb7i88++yzg5wQLZz3pEgPx3NxcLFiwAFVVVWNKpo8l0u3v70ddXR3sdjvmz5/vtwcoMPZIVxAENDQ0oLe3F/n5+QFF8/5CXkhzOBxoaGhAfE8P0qKjoVapIKhUICvnZ84EFRkJ6vhxCImJYkeazAGLee89qHt7QavVYN55B5TDIXrSctzpQpfsdUnFGgCUBQUQAMSXlEiWhPZ9+6A0m0FzHASdDti+HbbISEQdPixu5ePiQHV2iraIBoMY9TIMuLlzQbe3g8/OltILAKB87TVxirDZDFqnEw13FAqk5OaCiouTplYIKpXY0eZSWwxLeQDgfvjDYe+l81e/8hrpysG8+SbotjbwJSViLnqMUCqVUtohw3Wz4HleatMdHBxEa2srnE4nVCqV2zy7qKgov75HZ0rzBRBcL91w4qwmXYqioNFokJ+f7/bGh9pTd2hoCPX19VCpVJg5cyZOnDgREOECgZMu8fE1m81QqVR+u38RBPKeEGlaY2Mjurq6UFBQgPQtW8A8/bRoWi4bS83/6legX3wR6OqC8wc/cFOP0DQNhevCEARBnOhrNELIyADV3AwhO1uKcn1B6hSLiED63r1Q6nTieVMUlK4xNrWXXor4vj6kNDdDoVAgyjXfjGZZCEol+FmzwP7kJwDcI0YC589+Buatt0QnMZsNUCoR2dUFfuFCcail0wlERYHZswd8fr7Xdaqeegp0bS3ohgYx0p83T9IBe0a4GR9/DOXJk3AuXTraR+F1vSPBc/tP07RErASCIMDhcMBoNMJkMqGvrw8Wi8XtscRbl0SdZIgod801Z1SkCwR/NNHf/vY33Hvvvejr68OVV16J2bNnY/fu3UE7/llNugCQlpYWNk9dk8mE+vp68Dwvtc8KgjAu74XRIDfBSU1NRXR0NHJycgL6okkje/y4WEiRbGhoCMnJyW6FQDnZksfyPC9u011bdaXrd4Krk6zv0kthrq+HMz0d5r/+FTErVwIUBdu//jXiTSNi6VLRyWvmTERdeCEwNATnXXeJxuEqFYTkZFEdcNddKEtKAn7/e6iXLgXLcbBHRkJVUyOasUdHo+vii2GbNw8ajQYah8Nte00m7DL//CdAZqRFRoKNigKVnw/h6qtBnzgBAODLyiC4+vy9pRfohgaxaMcwoA8eFI29Xd9D0m2n3LEDET09YKqrwbz8MpzLloFubgZ9/DjYq68Ge+utI34+6gcfBADYN2zw+v/+5FxJnlitViNZ5v/AsizMZrPbxInc116DUqmEpr0dDMMgo7sb9NNPj3j88YBlWb/166Ey97nuuutw3XXXBf24BGc96XpqSsfqqTsS6cp1vsXFxW6J+/GO3hkJ/f39qK+vR2xsrGSCo9PpAo7k/SVdki6Ji4tDTEwM8n1EdYIHsSrffhsAwLkIQ6FQSBMn7HY7iu+/X3QD43lYNm0SPy+Ok86foijpD03TUG7bJo7GMRig2LNHfIzNBuW774oTg9vbIWRkoPGZZyDXFvDz54MGwJw6BdpiES0e//pXaCwWUB7ba7VaLUV0MTExEN57D8o334TS1apsVSqhMhqHqRSUmzeD+ewzUHq9eyHt4YdBHzwIimHAl5eLvr0+EHf8OGidDqAoqJ94AoIrH083N0O1bp3P3K8/GE+hi2EYxMXFuRnIqD79FE6nE0JXF5xOJ+x2O+qPHYNKpXJ7/6KiooKyzQ/ES9doNPrVPn+m4awnXU8ES/oFnC5W6fV6FBUVITk5OWh31pFI12AwoK6uDgzDDJtsTLb+gXzBR2sfNhqNqKurg0KhwKxZsxAREeF19DS5uRG1BiFJOViWRXNzMzTPPQdtYiJUv/+919ckpM3zvCi3EgTYHnhAPDeOA2bOhHL3blF5kJwsWibq9RCys8Fdey0EL62f9PHj4g9RURLxeXPBEgQByvvvFwtPq1ahp6cHVqsVitJS5F12GdRqNUyNjdCQvLaLXAn5ChqNOFH4d79zI17bBx+4ja1Rbt8OuqkJwOlCnHPpUrDvvAPBaARiYsQxSgwjjvQZBeoHHwTtGmapfvBB2DdsQISrWEny1MFWFzgefRSAmF5QAei4/HJUlpeDZVkpz97f3w+LxSK913K/hECDn0CmRgwNDZ11vgvAJCBdf7vSRoOcBJ1OJ5qamqDT6UJWrPJGuiSittvt0Gq1Xi3ryPMC+TL7anYg0ajFYoFWq5W+wCSClYMQJNnSkT+KN98E1dYGQRCg37IFJyoqMGvnTsQZjaBMJvAu717WZeAtXxMAqJ99FoqPPhLPjWHAZ2RAiIkBn58PduZMsPPnQ/XPfwIAnDffLNooUhS4O+4AXBNoCbgrrgD93Xeg+vvBFxf7HCET8cADUHz/PVQApv/hD7Bv2iQe3+kEtXs3HEYjIlpaMNjRAcuKFRj45S9FGZvdDpVKFXCXmCf6L7gACbW14jSMRYtAd3aKRjgBRriqjRul2XeqjRvhWLUq5KOjyNQLlUqFxMREt+5OItkkeeLGxkZwHCcV90hUPNJAzMluYA5MAtL1xHimRwiCIBWPpk2bFlCxKtAtv5x0HQ4HTp06haGhIRQXF7vl2TzhT1eaJzx1t2TmWl9fHwoLC5Gamuq2dvnPcrIl/+d5nhazGQODg2AKClBZWYmIL7+U7BdHA/3vf4NubgYcDjDvvw/+2mvF/2hshLBoEYRf/hI4eBACAMc990jr4O12pD35JCIjI8Ft3gzl3/8Oeu9e0VUsPh789OkBvUeAeMNWRUYiKjIS/QyDlORksFlZiMzKEvP5S5eKec6uLiiVSgwuXQpNf7/kgAW4zwpz/vKXkl+w3EGta/FiFHV0iO9vcTGoujpphLnqmWcAAI6VK6XiFYk27Rs2uOV0I267TZx8AYD+9lsACHgX5C/IGnDokM/HKBQKxMbGumnFiabYaDR6HYgpL9rRND3pDcyBSUi6Y4l05coAiqICGqcOeJ8C4c9zWJaVJG/5+fmYPn36qMQ9nvZh+STjnJycUW8qHMeNSLYmkwl1paXIaG9HamEhFK6mAvaRR6BybXs9I1xPUGaz2MDgMh8n4BctgpCaCmbPHjjffBMAQMpfQ0NDGNy+HWlGIyIEAeyHH4KHK/3B8+Bd+WKyds9ztG/aBPWKFdLPcpARN5a2NvA5OWB//WvEAoiNjZXcyJCZKXY+RkZicHAQbW1t4tQMlwyLRHSxb74JUJTbEEmyg3AuWSKuracH7JVXgk9Lg+qZZySvh8irrpJyvaonn3QjXgL2//0/KZ9ORsmHunkhUMg1xd4GYppMJrS1tcHsunnYbDZ0dXUhNjYWMTExI2qKpyLdCcJ4PHXlxjBEGeCreDQSSHTtL+nKSZ5hGJ+twt4wlkiXpmnodDocP34cycnJI3bLkciWZVl0vPIKoqKiwFx1FaJd5i7A6Vy3xWJBcXEx4srLR3x95k9/AgCwZMaZDHx5ubRFhkIBes8e8BdcII5T7+kR1797N/grrhD9gF0FzfK2NqgFATAaoaiuhn3jRtHS8dgxcAsWgF20CIIs9wy4F+wIfE1gaL35ZqRVVp4+hzffBP3dd+BnzhSbKCAqZ+R5YiLDUr7yCkx2O2wtLaBoGnHvvANuxgzYnngCkZGRYi68uVk8/7w8AAD34x9D4crX+gu50oH87It0lS+9JJ6na0DpRMOXyfnBgwehVqvR39+PlpYWr0VP8h4aDIapSPdMgD+euvI5aMQYRqVSob+/f0yRQiDyLzLsMiUlBdHR0QF7gQYa6er1euhcblYjzXvzLJLNmzcPju5ucaJGSwtMJpMb6ZLpxJ7krXjvPQAQhzYeOADlXXcBru0m86c/DSNe54YNoI8cEacJZ2eLv5Q1Eyg+/xx8TAxOabXo6elBQUEBUlJSwNTVAS7tJD93rugXsWQJhCVLxHHtrufLC3aCICDC5T/hnDEDAEA3NYECoHjtNTFP7AHmzTfBfPQRqKEhUP39UHz7LYSUFHDl5ZIPLuAuw1ISNYLDAWrfPghGIxy1tbCvWoXujAxoDx3CUHo6lCoVaKtVSk04Vq48nV7Ytm1YesEbPCVm3r6/ypdeEhtAXD+PlXiDZT/pCzRNQ6FQSI0d5DXtdrtUtCNFzz/+8Y8YHBxEdnY2qqurcd555wVsleoNq1atwkcffQSVSoXCwkK89tprQY+mz3rSDXRO2sDAAOrr673OQSNRciiMzMnrxsTESPKvvr6+gF4H8D/SJWbpTqcTSUlJyMzMHJFwPYtk6i++QITFglgA6f/8J3R5eTiRk4PExERERUXBaDTi8OHDACBFIRqNBglj6Fiyf/YZlHfcIeU14XSCrqsDVVsLR38/BhwOJO3YgdzHHjutGb75Zihc7b/czTf7PDZ5vEKhgHLlSlFzKwhgamrgcBmEC4IAgeelHRJFUcj88EMwR48CAYxuIpDmxL38Mphjx0AlJCACQER3N+LT02G1WhHV2QmOZdGfng6dTofWgwcRGRkJzU9/KkZ0djswAtn6QijTCxORuqAoSmqOkdc6ioqKsG7dOgiCgFdeeQUNDQ34/PPPx13wvuyyy7B+/XowDIOHHnoI69evx9NB1iWf9aTrCV85XTKah2EYzJgxw2u/NiHPYI7skcuxAp2/FuhrAae9bQcHB6WiHJmH5gl/imSOI0cwODQEgWEwPzYWtIdxDs/zUhTCbd6MDrsdbddcg5y2NkQuWgT8938j6a23xBy2l/SCdJxrrwX1xhug2tuB+HgxuunrQ6TDgeTkZNCJiaBckZ9ASN21S2DWrx81d+wGioIwYwbwi1+Afu01AAD3s59BKQhi9X3jRqR9+y0U9fVwXnIJHIsWIeL550EbjbDfcw+EuDhQMpcwX3D+4hfiSKLVq6H46isICQkQAJimT4fGaoUCcvOQWwAAIABJREFUQIJLcpYmCLBarTAajV7zxOTGRrbWPt9HL8TovOeeoKQXQt0CHEgknZiYCKVSiSVLluBHrrFLwcDll18u/Tx//ny8//77QTs2waQjXc+oUz4tori4eMQcUDA1vlarFfX19WJzQHGxzy3KWBsdPEGmRHR0dCA3NxclJSXScT2f4w/Zmi68EH1vvIEoQUDCwoVQq9UQGAaer0zTtFSxZtLTAQAZ8+bBMnMmjEYjDH19sHV0gOM4DHz7LTQaDWJjY6HRaNxubhxxdvrXv9CZloa2RYugvesuRP35z6L+VqYbpl0j4gXZRAF/4HzmGShXrpR+ZlzWjuxdd0HgOMlsvrK7G5E0DcFqhfLkSTC7doHS6wEAzMcfg8vNBbdgARwOB2ialkjOZxTI86JFpCCAOXwY9vnzYfOQh1EUhaioKERFRXnNE5tMJvT29op6YoXCjYhJ5Z88x9s6gpHLDcfUiEAi6VAX0l599VXc5DLuDybOetL1pdMlBuZmsxnFxcVep0V4IhhG5g6HA42NjRgcHBy1oWKsqgd51ErahBsbG5Genu5VeUFI17OTzBvZOp1OyYu46IYbkJiYCHrPHggA+Msu87ku5qWXQLnyhqpt26AaGEA8APbRR4H16yEIApKsVhgMBgwODoLduBEsy6Lvttug0WiQsWsXrFYrWLsdCUNDyDxxAmhrg5CdDZ6iRA8EoxFCVhZ4l+E5iXh9RblKV6uuUzYGyOnKmTLbtknrtTzzDL794Q+RnZ2NyspKKHt7Ibhy08IVV4gyKRfp0iYT+NRUsPfdB8b1fpLPw1eHnRwC/O9iHKldlxCxvPIfHR0Nh8OBwcFBaDSagD0MRkM4zG7OFFvHdevWgWEY3DpKW/ZYcNaTrieIQPubb75BYWEhUlJSAvLUHYuPAsMwEtl2dXUhLy/PLdL0hbGQrjxqHRwcRF1dHTQaDebOneszLSJNBB6hk4woKkikLJ90MRLZen293bsBV3sm40oJsI8+KkVy6enpYLKzIQBImj4dHR0d6OrqgkKhQCQAg9EIY08P1ByHiMhIRADgLr8cdFUVhOnT4XRNQWaeey6gdXmCEBQXGYmKigpJnsRdfbX0GO7qq8FdfTUiLrhAfM7ateJre3xmnjc18hlxHAcqPx/q774DN2MGBn/zG/Q1NyOwGN0dDMN4rfybzWYMDg66NSZERka6RcUqlWrMec9w2DoGcvyxqhdGs3V8/fXX8Y9//CMoOWJvOOtJl7wpRPDf29sLhUKB+fPnh9RpjIDneUkhkJ+fH5DGdyyaW4VCAaPRiCNHjgCAz/w0gSAIUCqVaGxshMFgkFIBxMqPKDkaGxuRmpqKefPmjenCYu+5B4wrbyjEx4OuqfH5WOYPfwB16hRYpxP6NWsQxzDQPPYYIiIioNixAyzLov8nP4HOaITBYEDMrl2gaRqJyclQG41gTSYk3nKLOKXBVelmXfpaAuXy5VLDgHL5cvdo1+lEzUUXQdPSgtS8PCjvu2/YGuXECwC26mrQu3ZJ8jZ61y7wP/6x9P/ygh0Bz/NgXnkFVFsbuFmzwE6bJnk+k7oDTdNeb4KBgqZpiVS1MrNzkifW6/Vob28fU55Yek/OMC9do9EYdMnYrl27sGHDBuzduzcoaghvOOtJVxAENDc3o6OjAzk5OeP21A1E49vb24tTp04hMjISmZmZyHNpLv1FoKTrcDjQ2dkJk8mE888/f8SUiTxvm5KSgvj4eJhMJhgMBvT29kq98g6HA1FRUSgpKUF8fPy47uysLG8oj3A9wXEcTIODcDidSM7LQ0REBFiXsoJbuhQUgGTg9JZ61iw4nU4YjUawGzfC/vjjsLS0IMJkAtraoNfrwRkMiImJGZG8BEFAR0cH2trakJeXh/QnnvB6voyri4yVdZGNBSTfS9G0ZCwuCIK0CyJpHnmKQhAEKBQKt3z8WOEtTwyIkj9veWK5CkWeJyY400g3FC3Py5Ytg91ux2Wu3d38+fPxkiuYCBbOetIlHS+eEeZYPXX9mR5BtvXR0dEoLy+HyWRCf39/wGv3l3Q5jkNzczO6u7uRnJwMjUbjk3B9FclIbjApKUnyXLDZbMjOzgbLsmhra0Ntba1UpCHFLm8Xnz/wRrZSCuOCCzDLbEZybCy4e++FP3sLpVKJxMREMK4ttSIqChTDQKBpCCwL9ve/x6Gf/hSAKGGLXbkS2evWQalSgd28GUNDQ6irq5N02d4ubkK2JNfLbN/uRrz8j38Metcu6We357oGWLIeAywtt9+OoQ0b4KRpJD70kFf1Cs/zUD77rNgye//9w9ITpDg2asHOTwSaJyZRscPhOOMMzIONBpdJfShx1pMuAGRkZLhV50PlqWsymVBXVweKoty29TabLSSeuiQya2lpQVZWFhYsWCDl7Lw9drQiGXEA0+l0Pqcjk4jSaDSiubkZZrNZ2rqS1MRYiFin06GhoQEpKSliCmPBAnieOePS3bIjjL5mXd4Dir/8BWBZIDcX8a4tZuIFF4hRtEvCdvK++2A0GmH54gvQNI3MzEwkJSWN64L1JFtfkEfV+b/4BXL//GfgxRel9RMwrsIemUgR9cILYF0KC3Lz9IyGfRXsAu1UdFvHCHlio9GIvr4+qXlIr9cHLU8sRyCRLvkMQ5FzDTUmBel6Yjyeut5I0Gazob6+HlardZifLnDaRyFQjES6ZJR6QkIC5s2bJ52Lt+cQI3VfRTI5AWRnZ2PevHk+SZNElPJImkRBRqMRLa7uNDJlgETEGo3G6zHNZrNkUzl79myvDRqEbKUIc+vWEYmXee45CAsXAseOQUhKcosuFQoF4uLioNFowLIsDAYDysrKpIYOkhIiQzDlUT1++Uswzz0HqrkZ/CWX+J1eYDZuBF1bK/08cPfdqK2tRVxcHH745Zegv/pK8lRgNmwYRry+II9uCUYq2BEJG8dxQc0TE8/a1tZWKJVKxMXFec0Ty9MT/uaJ5QiEdE0m01nppQtMUtINlt6WyKf6+/tRVFTkUwkx1nln3p5nNBpRW1sLpVLpNkpd/hxyofmjt+3v70dDQwMSExN9bqtHg7coiOM4iYjb29thMpkAQCLiqKgo9Pb2wmQyjahTHg+E888HX1Q07PfE/D01NRWVlZXSllWj0SAzM1N8rmwIZn9/vzjF+Y03kHz8OCL7+sDZ7bDffjvUarXf5MHzvFiUrK9HaWmpmGMeZbtMoloS8ZJ/+4Kvgp1Op8OpU6eQnZ3tRsTy542XiIkawleemOwwvOWJiaPYSK9PRsv7A71e73Xy9dmASUG6wfLUJaTLcRxaW1vR2dk5TD410vMChZx05dF0SUmJz6osiWT8cQAjHXjnn3++W7tzMKBQKBAfH+9GpoSIW1tb0dDQAKVSCaVSiY6ODkk5odFohuXtSFTrT3oBEJUKRC4mVy2QcfQ0TfuMqgnk5ubprqYO5auvQlAowCUnwxYdjZqTJ2Gz2aBWq90iYs8ozrlyJWyPP47BoSFQv/41ytPTT6tqXFEt43IH8xXljka2vmC321HrirLnzJkjnbM3XfZ4C3Yj2UbKawbSOcmMzjs6OqQbszxPLNcTB5LTPVttHYFJQrqeGA8Jms1mVFVVITMz02/513hez+FwoL6+Hn19fSNG08DpbiNSFIqNjUVcXJwbuRBvXtIUEs4v5tDQEBoaGpCUlISysjIwDCO1CRsMBnR1daGurg48z0sXHCEyhmFGJVtf4GTdZP42wng9zpIlYHbsAA0geulSzJkzB4BIbAaDAUajEd3d3bBYLOLcMFc+s7e3F7G33IKioiKfKS1/Uwr+QhAESVdNmnDk8BURk7+9pScAuEXD3lIbgeyUGIYZdmP2zBPL9cQ2mw0qlQrR0dGj5onP1qkRwCQh3fFGukSr2tDQAKfTiYULFwY8mSHQ4gzP81JRrKioaERvW3nEwjAMKisrxRZbgwHdLicwpVIpGUaT5oxwmZPII0zPqFreJkxALjyDwYCenh40NDSA4zgpAiJEPNJnwP7616Jsz+XalpWVhcrKynGdM/+jH0lKCl7Wz69Wq5GSkoKUlBTpd1arFXV1dejt7UV0dDT0ej2++eYbKac51oKjPzAYDKipqUFCQoJb+mQ0jESm8ojYV8GOZdmg54mB03riEydOwGKxoKamBnbXlA7POWzkWj9bvXSBSUK6nghEb0uixsjISMyZMwdHjx4NuAAXCOT63qioKGRlZWHatGkjPt6zSEZcl1JSUtzagOPj45GQkID+/n60t7dDpVJJhBcbGzvimJSxgIw1IhMv/G3J9Hbh8TwPi8UCg8Eg2W6yLDuMiEnXmMlkQm1tLSI9uslCDfJ+Nzc3Izc3F+eff75bg45nwZGiKDcijomJGbPsipjeG41GlJWVjdgUEwj8KdiRmzzRdpMccTAKdkRPrFQqUVBQIH2W8jwxGROvUCjw5ZdforW1FTExMbDZbH7ngUfDY489hp07d4KmaaSmpuL111+XagDBBDVKhBZaA80ggeM4t+19b28v9Ho9iouLfT6H5DwFQUBxcbFEAF9//TV+4OGk5Q/8eZ5er0dtbS2ioqJQVFQEs9mMvr4+TPcyWsafItnQ0BDq6+uh0WjcvqwEZFtMtsZWqzUoRCxXQ0ybNg2ZmZkhke6QYpf8HBwOh0QEubm5SE9PDxvhEqInn58/N2ciYSPrN5lMUnpFrvwYadtObtSNjY2YNm0acv/4R1AAnOvWBfHsfK+/qakJAwMDKCkpQUxMzLDUhBzjKdgdPnwYs2fPHvGmxLIsDhw4gB07dqC9vV0Krv7973+P22uC1B0A4IUXXsCJEyfG0xjh84KYFJFuIJ66pDHAbDZDq9WOyTAjUJDtN8uyKC0tlQjem77Xr5lkLq9cQRBQVlbm0y7S27ZYTsRdXV0BE/HAwAAaGhqkra363XcBANwtt4ztzRkBpNgVt307BAAtP/0pWltbkZWVhYiICJhMJhw/fhwOh0PyGCDnEKg950hgWVYyASopKfFZNWdco39Y1ygg4LSETZ5b95VeiYqKciNilUoFq9WKmpoaKJVKVFRUIHrtWtDffw8AUK5eHVLiHRwcRG1tLTIyMjDXZRRPzkl+LsEq2PnjMsYwDBYuXIivv/4aS5YswY033gin0xkUcx/550pGd4UCk4J0PeEtpyuf8FtYWIgZM2Z4fVPJEMdA79TenufN21YOuXrBH7J1cwArKhpTwWisRCwIAurr6wEAM2fODFlfuieY55+H88QJ6IeGEK3T4cLERPAyUgNODz80GAwYGhpCW1sb7HY7IiIiJBIjRDws/++yWPRGXvIIMycnB8XFxT4vRGbTJtCu94fZtMmNeD3hK69JCkw6nQ5NTU0wm83geR6pqalITU0N+eQGAqfTifr6ethsNsyaNWtE5ct4Cnbe0hP+Ep08pxvMdODq1auxY8cOxMXF4YsvvgjaceWYFKQ7UqQr95n1Z8IveW6gW1ZCoKQzqKWlBZ2dnSM6jpHnkJytr06ykRzAgoGRiJikMOx2O6Kjo5GUlASj0QiKohDzt7+JpuMAFG+9FfRo1263Q9fdDZXBgITERChcNy3PTa3yxRehBBC5bJmbF638ZtLR0SHJvwgRp2/aBKq2FhSGk6/ZbEZtbS0iIiLCkjOmKErSskZGRmJoaAjZ2dlIT0+HyWQ6fTO56iqU9PZCpVLBsnw5NFZrUHP1ZJxUXl4e0mXSt0AwUsHOV2OHvIjnT3rCYDCMqZA2mq3junXrsG7dOqxfvx5btmzBWperXDAxKUjXE0qlEg6HAx0dHWhubkZGRkbA8q9ALzJSvNPpdGhsbPTrNWmahslkQk9Pj9dtfbAcwMYClUolnU9eXh6ysrLgcDgkEuvs7ER8QwOi+/uhVqvBxMaCDhIBkBtld3c3ClasQMJbb4ECwN5//7DHMlu2QPHBB26/Y5ctcxvzIp9CS8aBGwwGDAwMQD00JI7+HhyESqWC3TUmfGhoaES9tCfYFSu8phcCgTzClM/8io6OdruZ0Dt3wuFwoMdoRGdnp7QzkRcc5ZV+f2Cz2STvjVDdZHxFxWazGTU1NUhOTvY7Ih6rl+5oto4Et956KxYvXjxFur7gSVQDAwNS4YIMnfQXSqVyTJpblmXxzTffuA269AVyl1coFCguLoZer0dHR4fblphhGPT29iIqKgpz5swJao5yNBBFR1xcHObOnStt34ZFxHPmQNixA1arFV2XXgpDTY2ktRxrsc5bN9lIz6T37pWmOtB794L/r/8a8fhy5Qe2b4dy9WpwHIeeBx5AU3c3+g4dAsMwiIqKQk9PDywWi5sVpj/wt7uMgEylbmlpQX5+PtLS0ny+FkVREJ56CkoAhbLfkxsi6Qgjlf7RzIvket/i4mK35oZQg9xce3p6MH36dOkGN1JETNDT0xN0yVh9fb1UfN+5c6fXAncwMCnUC4C4FdXr9airq4NarYbBYMCFrgkDgeDkyZNIS0vzO19KTHCMRiO0Wq3bJFNPjJa3FQQBer1ecgBTqVSScFxOYqHa6pIRQxzHQavVjnmem3xbbzAY/CJiq9WK2tpa0DSN4uJivzvomOefh2LnTgDiyB9v0fBosFgsqK2thUqlQnFxMVQqlTQmh5yDxWIBwzBuOWL5WHpAZl7j8lngtdpRiZdEeYEoIvyF3LzIYDBI5kVEOcEwDNra2hAXF4fCwsKw7aIAsd395MmTSEpKQn5+/qjpBHLN2Gw2bNq0CTt27MCJEyeC2gp8/fXXS9/B3NxcvPTSS8jKyhrr4Sa3ekEQBBw7dgxOpxPTp0+HRqPB119/PaZj+dtdZrfb0dDQAJPJBK1rPLgvMvSnSEYcwPr7+1FQUCCN+ZEXiQYHB9HS0iJ54MpJbDwXK5EF6XS6oEQ7oxXrOjs7pfxqTEwMrFYrLBYLtFptwK8tJ9lACZdYZup0umFKFpVKhaSkJLf1EBIzGAxobGwcFk2mORxQ+fk5kNfu7+8PKI0RCLyZF3EcB4PBgKamJhgMBqhUKgwODqKmpsZNSxzsUT8EPM9L3YNyJc9ooGka33zzDZYvX46f/OQnaGpqCrqe/q9//WtQj+cLkybSHRoacouOvv76ayxYsCDg/GJTUxPUarVPUbTc27awsFDaCjY0NECj0biZgMglNL6KZIIgoLOzE62trcjOzkZWVtaod31v+lWn04no6GiJhEfr6CLHIUL/rKwsZGdnh62LjWh9m5qapIjR0+fA39QE8/zzAAIjXdKAkZmZOa7zJk5m8mgy7513oFKrYb3vPqmbSn58kkLJyMhATk5OWMeaExlYeno6pk2bJhV+5Vpio9EInucD6hD0B3q9HjU1NUhLS5Ne2x/Y7XZs2LABX375JbZu3Yrzzz9/XOsIE3x+aScN6TqdTjex9sGDBzFnzpyAvyhtbW0QBGFYl5int63nl8aTrD07ybwRh9wBLC8vb1xfak8iNhgMUmutnIhJBENSMRqNBoWFhSHtwvOEvJusqKjIbYfgLTUxFiL2BZLGUCgU0Gq1IcmVE6MXcg7ECjMyMlKKjondZLjgdDrR0NAAq9WK0tLSUdM38g5BckNhWdZNS+xvqovjOEnqWFpaGlDa6ujRo1i+fDmWLFmCVatWhfV7Ok6ce6R79OhRlJSUBPzFJjrVgoIC6XfElyExMREFBQVeP3hC1jk5OaOPN5c5gBUVFQXdAYxAfuGQPyzLgud5UBSFgoICpKamhi2XR3TGBoMBWq3W7y21P0SsefZZUBQFp8docwKe59Hc3Iy+vr5xmeKMBYIgoK2tDa2trUhMTIQgCMOsMIluNxSfBWk7H48MDHC3wyRE7HA4RtRDDw0NoaamBpmZmcjJyfH7te12O5566ins378fW7duxcyZM8e05gnE5CddYslIcPz4ceTm5gacaO/r68Pg4CC0Wi0MBgPq6uqgVCqh1WpHJMeOjg5YLBbkueakeSNb4gBG8sDhdADjOA4tLS3o6elBZmamNOCSbCXlBDYefwBvkKdQcnNzkZGRMW5ZmZyI4557DupTp0ArFOBc0Sv/+ONSREymVsi31OGC0WhETU0N4uLiUFBQ4JYrlU+5IBGxvEV4vPlVuQxMq9WGpABLag7kHIxGI2w2m6QC4nke06dPR0JCgt+f+eHDh3H//ffjhhtuwMqVK0OWXw4xzj3SrampQUpKSsCFmcHBQbS1tUk5xtHIkeRsyZaZtKTKi1wKhUKSxuTn5yM1NTVsY0YEQUBPTw+ampqkaMObaF2+HTYajQAwjIjHQlbEbyI+Pn4Y6QQLynXrQNXUgOc4WPLyxCLnbbfBarVKLaK5ublISkoKuumPL3Ach1OnTkGv10vFXX8wUn5VHk2O5tUwUTIwQEyb1dbWIjExESqVShyZJLPDlGuJ5d8pm82G9evX48CBA9i6dStmzJgR1nUHGece6XorbPlzjJqaGvT09OC8884b1dvWW5GM3Pn1ej0MBgN0Oh2sViuioqKQnp6OuLi4US+aYIFE6tHR0SgsLAwo0pEbtRACCGRWGlF32O12aLXaoDli+YLS1UnmXL1a6gjs6elBXl4eFApFyHLE3kC28zk5OcjKygpKswhpESafhdyrQV7oMplMqKmpQWxsbNhlYKS5w263o7S0dJj7l9PpdLuZEAnb22+/jejoaHzyySe4/fbb8eCDD56t0a0ck590PZ3GWlpaoFAokJ2dPepzeZ5HW1sb2tvbkZGRgcHBQVRUVPh8vD9FMrkDWH5+PliWlYiYXDTybWQw83mE8EikHqxZUmQyhGeBSB7VR0ZGor29HZ2dnSgsLBzxxhUKDAwMoK6uDmlpacjNzfV6UwhVsY5s52maDlmRjoB4NcgLXRaLBYIgID09HSkpKW5WmKEGqXsEmjc2Go149NFHceLECaSnp6O9vR35+fl412WkdBbj3CPdjo4OOJ1OKcfqDXJv25SUFOTn50MQBBw9ehTz5s3z+vhAHMCKi4t9VmpJ9CInYmB8W3oyZqinpwcFBQVhITziIUt8cPV6PVQqFZKTkxEfHx9wN9dYYbPZUFdXB0EQRs2/e8N4iJjctLu6uiZkOz80NITa2lqkpqYiOTnZLU/sdDrdBnAGu7nG4XBI00BKSkoCutFUV1dj5cqVuOWWW7B8+fIxje05gzH5SZfneTdnsdE8dUmrK+kEIlshQRBw4MABN29cfx3AiKH3WB3AyJaeELG3SNKzC4qsr7e3F01NTRNSLCJTFCiKQnFxMRiG8dnNJY+Ig0HEcsLzNrZmPPCHiO12O+rq6pCUlCSlMsIFuQxs+vTpXpU6ZCqDPCL2tMLUaDRjMgIn5jgFBQUBpfGsViueeOIJHDlyBNu2bUNJSUnAr30W4Nwj3YGBAfT09KC0tNTtccTblrS6ett6E0Nyf8jW0wEsGJV5OUgkSYiY6DzJfDSFQoGWlhZERUWhsLAwrB4N8o6u0WRY3tpqx2uoPjg4iLq6OiQnJ4eN8OTua93d3XA6ndBoNFJUH4ocsTeQHdpYvnPyLkfymRDfD09PYm/HtdvtqKmpgUKhQElJSUDa2QMHDmDVqlW4/fbbcd99902GiNYXJncbMDD6nDQi1xoaGvKr3VTuBeqrk4w4gKWkpITMAczb+HOn04n+/n60tLRI89EoikJbW9uoF0wwIPeZ9Xc2mbe2WrlrGdFHE+vFkc7Dbrejvr4eTqfTzY0rHCB+GDqdDkVFRUhPTx/mvhbKYp08bzxWNzCKohAZGYnIyMiArDA1Gg0GBwfR2trq1R96JFgsFvzud7/DsWPH8M4770Cr1Qa87smCSRPpCoIAh8Mh/dtqteLkyZOYNWuWNE49Pz9/1KiApBdSUlKkaNLzTm4wGFBfX4+IiAgUFhYGbUaTPyDOTF1dXVJzA0VRbheMXq93cywj5xGMXB6RxkVEREjmMMGEty09OY+YmBiYzWb09PSgsLDQzbIxHCBDE8m5jxThBbtYRzoi29vbg55GGQnkPAYGBtDd3S1puuPi4qRzGSlVRK6nVatW4Y477sD//d//TeboVo7Jn17wJF2n04mqqirQNI2MjAzk5uaO+GHLUwl2ux1DQ0PD2mmjoqKg1+ulYk2wVAH+QB5Zk9710c6HbCHJH7lRDrlo/N0ajrWbbLwgEVhXVxfa2tqkIYpyLXSwbii+IO9mKykpGbOloC8ilutvvRGx2WzGyZMnodFoUFRUFFbSkpO9VqtFYmIi7Ha7m3yNaHDl+lsyNHLt2rU4ceIEtm3bhqKiorCt+wzA5CddQPxSA6dlQxaLBRdeeKFf3rYj5W1JakKn0yE6OlpKW5ALJS4uLmTjtoHT9pEqlcqt6BcoRvNnIEQsv6hD0U0WCBwOh6T9LCkpQXR09Kg3FHLxB4OIyXcpVAXKkYg4JiYGer0eQ0NDbn6z4YLFYsHJkycRExMzKtnLc/ZVVVVYt24dzGYzysrK8D//8z+45JJLkJOTE8bVTzjODdIdGBhAbW0tKIqCVqvFt99+63NCrz9kO5IDmDfNKilwEQIbb4U+HG3D8iGJJHIhW0iVSoW+vj4kJiaG3RSHdFW1t7e7pVFGejyp0pM/RC41FgtMIoUidqGh8sfwBhLZt7S0SPl6fyLiYIF4RXR2dmL69OkBRfYmkwm//e1vUVdXhw0bNkCv1+PIkSOYPXs2Lr300pCs9wzF5Cddoq/NzMyUik7exqL7Y7cIjM0BjHTcyCv08sJQXFycX+oCuSJitEkCoQCZQEs66Ww2m1s3GonsQ7Um0jqckJCAgoKCMW+nfUX2nkQs736S32j9Iftgg2VZaVp1aWmpVCQca2oiUJjNZpw4cQIJCQnIz8/3+70XBAFfffUVHn74Ydx111245557wipbPAMx+UkXGO405ump608nGXEAI6N0xhvhkIkW5GKx2+0jRl/EnCUlJSXsuk+iefXWTSZvgiC+sfLIPhhNEA6HQ9KdlpSUhKR1WN7J5dlfvA+lAAAZyElEQVQdqFar0d/fj7i4OElvHE4EKgOTG82Ml4hJ63Rvby9KS0sDMooyGo14/PHH0djYiO3bt4/YkHQO4dwkXeKpyzCM33lbk8mE4uLioM9fIiDbYDkRcxyHiIgIKTKWRzjhAonsA9G8eovsx6K9JcWatra2CYnsnU4n6urqMDg4CI1GA7vdPsztK1S2i8Bp3StN0ygpKRlXLnosRExG55DP3t8IVRAE7N27F4888gj+93//F3fddVdYotu2tjYsXboUPT09oCgKd911F5YvXx7y1w0Q5ybpHjlyxK0H3hvZchyHtrY2dHd3Iy8vb0IueKIfTkpKkkbCkJHcxCAnJiYmJOsi3WQAxtQ+64mRJF9y7S2BwWBAbW2tV+vDcIDsLMjkDPIe+3Je8yTi8ZCM/GYTqO41EPgiYo1GA6vVCqvVirKysoDUOEajEb/5zW/Q2tqKbdu2ITc3NyRr94auri50dXWhvLwcRqMRFRUV+PDDD1FWVha2NfiBc4N0idMYydk2NTWhu7tb0hXKLxRiedjc3DwhrbM8z0tSHG/bSXmhTq/Xu23nybmMp1An7yYrKioKmV8AkXzJI3tifO10OsFxHEpLS0O2s/AF4tUAiDcbfxQhns5rJpNJujkG6pchl4EVFhaG/WbT19cnTe+gadrv1IQgCPjyyy/x6KOPYtmyZfj5z38+4bnba665BsuWLcNll102oevwwLlBuk6nUzJOJnlbksPT6/XQ6/VSdd7pdCImJgaFhYWIi4sLa3RLtvKkX9/fC06+ndfr9W4dXISIRyvUybvJfPnrhhKkUNXc3IzExETQND1Muib3IQ7F65O8dTCaDHw5r/mywJRrfidCBkZ8fg0Gw7DROSOlJr799lsUFBTg1VdfRWdnJ7Zu3TpspNVEoLm5GRdddBG+++67gPLQYcC5QbqrVq1CTEwM5s6di4qKCmg0mmFjvuvr68GyrNS+qdfr3VQG/pLXWGA2m1FfXx/wmPGRINer6vX6YQ0Q8uo80fuq1eqQdJONBqPRiNraWulmJy8gym+OntI1f7fzysceAwA4f/97r/9vMBhQU1ODxMTEgCrzgcJX0VGlUsFgMCA1NTXsTQ7A6aGUnqmUkUC+X7/5zW9w4MABOBwOlJWV4eKLL8bDDz8chlX7hslkwn/9139h9erVWLJkyYSuxQvODdKtra1FVVUVqqurceTIETgcDsycOROlpaU4cuQIFixYgBtvvNGrMYt8Cywnr7i4OKlpYKxbQLkDWXFxsZuPQrDhTSbFsqyUciGt0OGMblmWlaKrsUxRkLuueea6SRSpfOwx0CdOiM8rK3MjXiLDMplMmD59esgN1T3Bsizq6upgMBiQmJgIm80Gs9k8zHktVBaY5PwtFotfQynl0Ov1ePTRR9Hb24uXXnoJOTk56OzsRFNTE374wx8Gfa3+wul04qqrrsIVV1yBFStWTNg6RsC5QbqesFgseOyxx/DOO+9g5syZ6OnpQVRUFCoqKjB37lzMnTvXp9E1IS95LlIQBOmC96cLTV4omTZtGjIzM8OaxiBb+ZaWFqSmpkKtVkuu/RRFhVx3Kx/xHqzz95Xrnv7664hpaQHDMKDOOw/s73/vlkqZiPcfOG3u7e31vak/SDttsCwwybj3nJycgM5fEATs2bMHjz/+OFasWIGlS5dOeO6WQBAE/OxnP0NiYiKef/75iV6OL5ybpCsIAl5++WXccsstUvvowMAADh06hKqqKhw6dAgtLS3Izs5GZWWlRMa+hujxPO9msyjvQiNETIoPAwMDqK+vl7ay4S6UkBHrsbGxXicYcxzndsF7Rl7ycxkLiDEO8SsOZTcbUXxErFkjyq/+53+gUCikgh3xyQgn4drtdtTW1gJAQObecscyg8EAq9U6JhkekcGRjrpAWseHhobwyCOPYGBgAC+99BKysrL8fm44sH//fixcuBDnnXeedCN48sknsXjx4glemRvOTdL1B6SwUV1djerqavznP/+B0WhEaWmpRMKzZs3y+aUl0QohYrPZDJZloVQqkZubi5SUlLDmTolXgc1mC7jBwNsFr1arpa28P8YyLMuisbERer0+7BOPgdMi/66uLmRmZko3F3lRKJR5e3lHW1FREVJSUsZ9zJG60bxZYJImi0A1z4IgYPfu3VizZg1WrlyJ2267LWzR7Z133ol//OMfSE1NxXfffReW1wwxpkg3EDidThw/flwi4mPHjoFhGJSXl6O8vBxz585FcXGxWyGEZVk0NTWhv78feXl5oChKImKWZREdHe0mWwt2EUXeOhys9lVfci+S65a30cq38sEayBgoyNga0s0nJwxPv1jSHUjcysj5jCciN5vNqKmpQXR0NIqKikK2u/F1LiqVCg6HAwzDoLS0NKAb7uDgIB5++GEYDAb86U9/QmZmZkjW7gv79u1DTEwMli5dOkW6wV/L2QdBEGA0GvGf//wH1dXVOHjwIBoaGpCamory8nKYzWYMDAxg7dq1XqvC3mRrJKcqLwiNlaRIKiMcExS85bpZlgXLsoiIiEBBQYEkBQsXyBRam83mc2yNN/gyyfGUro1GnvIW2omQgZHcOTHUpyhKukHKLTBjY4fPRxMEAZ988gnWrl2Lhx56CLfccsuE5W6bm5tx1VVXTZFu8NcyOSAIAj766CM8+OCDSEpKglqtRl9fH4qLiyXJWnl5uc+KNCkIESI2m81QKpVSbpjk7kZCsLvJAgXHcWhsbMTAwACys7PB87wk9yJa1WDcVHxBEATJjStY7cOe3gwGg2HElmC9Xo+amhqkpqb6LMqGEqSFmGEYaLXaYTI8+U3FaDRKO5X9+/cjISEBu3btAs/zePHFF5GRkRHWtXtiinRFTJHuCPjXv/6FnJwcafglx3E4efIkqqurcejQIRw5cgQcx+H888+X1BKlpaU+IyeiGyY5YrL9lRMxwzBSNxkh+XBPnyWG6qdOnUJ2drbX6N5zthu5qci38uMp1Mm38qG2nfTWEix3qiMtvOFuMiE3nEBaiMlO5ZlnnsGePXukvH1JSQnefPPNsBd85ZgiXRFTpDsOkEjj8OHDOHjwIKqrq1FTU4O4uDhUVFSgsrISc+fOdfPp9fZ8T5cyp9OJxMRE5ObmIjY2NqwXu8ViQW1tLZRKJbRabUBFQlKoI+cj92UgRDza8TiOQ1NTEwYGBlBSUhL2rTwgysDq6+ulIqnRaPR7cnMwQEZRRUZGBuyG1t/fj1WrVoFlWbz44otIS0uTutQmem7ZFOmKmCLdIEMQBOh0OqlId/DgQXR0dCAvL0+KhsvLy4e1JsunR2RkZEhkTLbycuIKhche7tWg1WqD0uDhawKE5yQLQipEc5qRkRH29mVgdBkYy7LDZHjB1N0SY/eOjg5pdE4gz/373/+OJ598EqtXr8ZNN90U9kLnaJgiXRFTpBsG8DyPhoYGiYQPHz4Mi8WCGTNmYObMmTh8+DAWLFiAW265xWtkJ7/YvbU1j3eGGEklZGZmIjs7O6Rk5y2nynEcOI6DQqGQzHkmwi9iLDIwXzI8T93taCCjc4hBTiDFUp1OhwceeAAURWHLli1hH+jpD26++WZ8+eWX0Ol0SEtLw9q1a/Hzn/98opc1HkyR7tkGu92OJ598Ei+//DJmzJiBvr4+REREYM6cOVJEXFBQ4JN8vHkyyCNIfwxlrFYramtroVAo3CwywwUS2bW1tUkpGHmhTt7IEaoW2lDIwDyje7l0zVNlIAgCWltb0d3dHfBQTEEQ8OGHH+Kpp57CY489hhtuuCGs0e2uXbuwfPlycByHX/ziFxPu1RBmTJHu2Yht27bhxhtvRHx8PARBwNDQEA4dOiQV6hobG5GVlYXy8nKpoy45OdnrheWrrdnbcE3SMNLb2xvwNjZYMBqNqKmpQXx8vNeRPd628vLOLTIaaawkI5eBjWcCsD8gaRb5Z+N0OqFWq2E2mxEfHx+wuXlvby8eeOABKJVK/OEPfwhKk0Yg4DgOWq0We/bskTo+33777TPN8zaUmCLdyQgyXqeqqgoHDx7EoUOHMDQ0hJKSEqlQN2vWLJ95RG9tzUR4n5SUhIKCgpBFkL4g72gLxBwHcFd/yAt18kYOf4iLyMC8NVmEA+Sm19PTg9TUVKnNmYwVGmmahSAI+OCDD7BhwwasWbMGS5YsmZDc7YEDB7BmzRrs3r0bALB+/XoAwCOPPBL2tUwQpkj3XAHLsvj+++8lb4lvvvkGFEVh9uzZUltzSUnJsIvVZrOhtrYWgiAgLS1NKtR5EldcXFzI5FmkfTVYHW3yQh0hY9L8ICdi8l4QNzTSBi73mg0XyOiclJSUYbpfMrlZbn8JABqNBvv27cO0adPwl7/8BdHR0di8eXPIJlH4g/fffx+7du3Cyy+/DAB44403UF1djS1btkzYmsIMn1/eiRPlhQCPPfYYdu7cCZqmkZqaitdffz3s7YwTDYZhMGvWLMyaNQt33303BEGAyWTC4cOHUV1djaeeegp1dXVISkpCRUUFZs2aherqapSUlODaa68dpvmVE9fAwACam5uD3tZMcscMw6CioiJoXhUURSEyMhKRkZFIS0uTzocU6rq7u1FfXw9BEMAwDMxmMzIyMjBnzpywe93yPI/GxkYMDg6irKzMawsvaTiRR/8cx0Gv1+O7777D9u3b4XQ6kZqaig0bNmDDhg3hPIUp+IlJFekaDAbJPf6FF17AiRMn8NJLL03wqs48kLbRV155BS+++CKmTZsGi8WCwsJCKS0xZ84cn3PZSMRFIsixtjXzPC8ViSaiyQMQUxInT56E0+lEUlISzGazm+Y2lDI8ApLOIGOjAnmd7u5urFixAhqNBs8//zySkpJgMpnQ0NCA2bNnh2S9/mAqvXCORLrycR1ms/mM0yGeKaAoCunp6XA6naiqqkJubi44jkNdXR2qqqqwc+dOrFmzBg6HA+edd55ExGVlZVAqlW4RF7H9k1tFNjY2unnDemtrJuY0ycnJqKysDHtkKZeBFRYWDpNRyQt1xACcFOrkLmXj+Y5xHCeZq5933nkBTYDmeR7vvvsunnvuOTzxxBP4yU9+Iq0lJiZmQgkXACorK1FfX4+mpiZkZWXhnXfewVtvvTWhazpTMKkiXQBYvXo1duzYgbi4OHzxxRdhr9pOJthsNhw9elTKD3///feIiYlxM4EfaaCnt7bmiIgIOBwOCIKA0tLSCekoI5rXQGVgcmcveZv2WFzKBgYGUFdXF9DoHILu7m4sX74ciYmJeO655yZEXeIPPv74Y9x///3gOA533nknVq9ePdFLCicmTyHt0ksvRXd397Dfr1u3Dtdcc4307/Xr18Nms2Ht2rXhXN6khiAI6O/vdzOBb21txbRp0ySTn4qKCq8m8MQroKmpCYmJiaAoSqrIy2Vr/k7THQuCLQPzlHrp9XqwLCspDEhHnacFaH19PaxWa8Cjc3iexzvvvIMXXngBTz75JK688soJ28299957WLNmDU6ePImDBw9i7ty5E7KOMxiTh3T9RWtrKxYvXhy0lsJVq1bho48+gkqlQmFhIV577bWwjw0/E0FG3ctN4E0mE8rKyqSImGEY7Nu3D5dffjmKi4vdokHPOWhGo1GaxkGIa7wjawAxb0rSGaGUgXkO2DQYDAAgDdXU6XTIz88PeHRQV1cXli9fjpSUFGzatCmkc/b8wcmTJ0HTNO6++24888wzU6Q7HOcG6dbX10uOX3/4wx+wd+9evP/++0E59qeffopLLrkEDMPgoYceAgA8/fTTQTn2ZIPD4cDx48exf/9+vP7662hvb0dZWZk0jaOyshJFRUU+iY/kUwlxkbbmQCZYEJC86UTKwGw2G06ePAmbzYaoqChYrVa3G8tIhTqe5/HWW29hy5YtWL9+PRYvXnxG1SoWLVo0RbrecW4U0h5++GHU1taCpmnk5uYGVblw+eWXSz/Pnz8/aGQ+GaFSqVBRUYGqqirccsstWL58OaxWq2QCv2bNGpw6dQppaWlu+WEy7YJhGCQmJrrlKolsbWhoCK2trVJbsze9LYFOp0NDQwOys7Oh1WonhKx8jc4hDQ96vR69vb1SoS4uLk6SsKnVatx3333IyMjAvn37pnZWkwSTKtINF66++mrcdNNNuO222yZ6KWctyKTk6upqKT/c398PrVYr5YfnzJnjMwKU622JbI20NUdHR6O/vx80TaO0tDTsnhGAGO3X1NSAoii/W3hJoW7//v3YuHEjWlpaUFpaiiuvvBJXX301zj///DCs/DT8qZ9MRbo+cW6kF8YLf75k69atw3/+8x988MEHZ9Q2bzKA4zicOHFCcls7evQoBEFwM4GfPn26T7UBy7Jobm5GR0cHoqKiwHGcNOE4GMbp/kAQBPT09KCpqcmrFG00tLe347777kNOTg42btyI/v5+HDx4EBkZGVi0aFFoFj0OTJGuT0yRbjDw+uuvY+vWrfj8888D0lT6g6lq8HAQkx65CXxtbS0SEhIkpURlZSWysrJw8uRJDAwMICkpyW3ku9PpdFMXyNuagzGMUg6bzYaamhrJ4D2Q4/I8jx07dmDr1q3YuHEjLrvssrPipj5Fuj4xRbrjxa5du7BixQrs3bs3JNrfqWqwfyCjgohaoqqqCt9//z1UKhX++7//GxdddBHKy8sRGxvrMy0h92PQ6/VuRjJE5hWIuoE0WrS1tY2ps66trQ333nsvCgoKsGHDBrcmnzMVf/vb33Dvvfeir68P8fHxmD17ttR9NgUAU6Q7fhQVFUnuW4BYTAtFi/FU5OA/BEHAVVddhXnz5uH666/HkSNHJBN4m82GGTNmSGmJmTNn+syrejOSoSjKL79eMjonKioqYL9dnufx+uuvY/v27Xj22Wfxox/9aMKi2ylJZNAxRbpnC6ZINzBYLBavqR673Y5vvvlG8h7+7rvvEBERgfLycomI8/PzfUa0pK1ZPliTqAsIGff19aGjowMlJSUB62ZbW1uxbNkyaLVabNiwwavBTTgxJYkMOs4NydiZDn+76abgP3zl1tVqNS644AJccMEFACCZwJPc8AcffCD5AhASrqioQFJSEiiKgkKhQEJCghuZEnWBTqfD999/D4qiEB8fD71eD0EQpGnNI4Hnebzyyit47bXX8Oyzz+KSSy45I3K3U5LI8GGKdMOIzz77LOyveY6PTJFAURQSEhJwxRVX4IorrgBw2uWsqqoKX331FTZt2iSZp3uawAOAQqGQRrGTvLHFYoHBYEBvby8aGhrA87yb25q8rbm5uRnLli1DWVkZ9u/fP+HRrS+8+uqruOmmmyZ6GZMWU+mFMwzBTC9MjUwJHE6nc5gJPE3TyMvLw/Hjx7Fy5UrccMMNPpUJ8rZmvV6PoaEhPPTQQ0hJSUFtbS2efPJJ3HzzzWGfRgFMSSLDjKmc7pmOUFSDpzxNxw9BEPDb3/4Wf/3rX3H55ZejtbUVdXV1SElJkbrpKisr3brN5GhqasJ9992H1NRUlJaW4ujRo7DZbPjkk08m4GxGRiglkecgpnK6Zzquu+46XHfddUE9ZkdHB3JycqR/Z2dno7q6OqivMdlBURQWLlyIxx9/XMrXEse0gwcPoqqqCtu2bUNvby+KiookIp41axbefvttvPHGG9i8eTMWLlx4RkeOu3btwoYNG7B3794pwg0xpkh3ClMYBZdddpnbvymKQmZmJq699lpce+21AMRUTm1tLaqrq/Hhhx/innvuwbx58/Dvf//7rCCxZcuWwW63S+caKknkFKZId1IjKysLbW1t0r/b29ulSQ/Bwp133ol//OMfSE1NDZqN5tkIhUKBsrIylJWV4Y477oAgCBMa2QY6L7ChoSGMqzu3MZXTncRgWRZarRaff/45srKyUFlZibfeegszZswI2mvs27cPMTExWLp06TlNumcapuYFTjimcrrnIhiGwZYtW3DFFVdII1OCSbgAcNFFF6G5uTmox5zC+DE1L/DMxRTpTnIsXrwYixcvnuhlTGEC4DkvcApnBqbSC1MYN5qbm3HVVVdNpRfCjKl5gWc0pnS6Uwgdpkj3zEaw5wVOwS/4JN3wt8VMYQoBoK2tDRdffDHKysowY8YMbN68eaKXdFagvr5e+nnnzp2YPn36BK5mCnJMRbpTGBduvvlmfPnll9DpdEhLS8PatWvx85//PGjH7+rqQldXF8rLy2E0GlFRUYEPP/xwqpV5FFx//fXD5gUGWy44hRExlV6YwuTANddcg2XLlg1rWJjCFM4wTKUXpnD2o7m5GUePHpXsGs81PPvss6AoCjqdbqKXMoVxYIp0p3BWwGQy4frrr8fzzz9/VoyzCTba2trw6aefYtq0aRO9lCmME1OkO4UzHk6nE9dffz1uvfVWLFmyJOjHt9lsmDdvHmbNmoUZM2b8/3bukKWhMArj+P9oGMqCRVBmMFgW5BaRhVsV9zmWDIatCDOs3+4HsFpkICadQavJYjSZtVgGPpYbBhuIuO1u8/nF94Zz0uFw3vdcOp3O2GP8VbPZJMsyLzksAC9H2EyTRKPRoFqt0mq1JhKjVCrR6/Uol8v0+33SNKVer1Or1SYS77e63S6VSoUkSYpOxcbgp4s0s0JFRAo8AM/AV37clnQzoXirwCNwLGlq/8GMiFtgY8SnM6ANHEr6iIhXYE+SB7tzykXXDIiIZeAJ2AHOJZ0WnBIAEbEL3AGf+dEW8AbsSxpeR7OZ56JrNiAi1oAr4ETSzK1wudOdf75IMxsg6R24B46KzsUWk4uu/XsRsZ53uETECnAAvBSb1WiStt3lzje/XjCDTeAin+suAZeSrgvOyRaUZ7pmZlP0DT58wkRt8APlAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Create a 3 X 20 matrix with random values.\n",
    "\n",
    "mu_vec1 = np.array([0,0,0])\n",
    "cov_mat1 = np.array([[1,0,0],[0,1,0],[0,0,1]]) # eye: 线性无关\n",
    "samples =  np.random.multivariate_normal(mu_vec1, cov_mat1, 500).T\n",
    "\n",
    "print(samples.shape)\n",
    "\n",
    "#Compute the mean vector\n",
    "\n",
    "mean_x = np.mean(samples[0, :])\n",
    "mean_y = np.mean(samples[1, :])\n",
    "mean_z = np.mean(samples[2, :])\n",
    "\n",
    "print(mean_x, mean_y, mean_z)\n",
    "\n",
    "mean_vector = np.array([[mean_x],[mean_y],[mean_z]]) # 3x1\n",
    "\n",
    "#Computation of scatter plot\n",
    "\n",
    "scatter_matrix = np.zeros((3, 3))\n",
    "for i in range(samples.shape[1]):\n",
    "    tt = samples[:, i].reshape(3,1) - mean_vector\n",
    "    # 协方差有正有负抵消, 方差一直累加\n",
    "    scatter_matrix += tt.dot(tt.T)\n",
    "print('Scatter Matrix:\\n', scatter_matrix)\n",
    "\n",
    "ax = plt.figure().add_subplot(111, projection = '3d') \n",
    "ax.scatter(samples[0], samples[1], samples[2], c='r', marker='+')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 99,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "40\n",
      "Scatter Matrix: [[ 369720. -123240.]\n",
      " [-123240.   41080.]]\n"
     ]
    }
   ],
   "source": [
    "arange = np.arange(0, 40)\n",
    "samples = np.array([arange * 3, arange * -1]) # 线性负相关\n",
    "mean_vector = np.array([0, 0])\n",
    "scatter_matrix = np.zeros((2, 2))\n",
    "for i in range(samples.shape[1]):\n",
    "    tt = samples[:, i].reshape(2, 1) - mean_vector\n",
    "    scatter_matrix += tt.dot(tt.T)\n",
    "print('Scatter Matrix:', scatter_matrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 123,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Scatter Matrix: [[47970. 15990.]\n",
      " [15990.  5330.]]\n"
     ]
    }
   ],
   "source": [
    "arange = np.arange(0, 40)\n",
    "samples = np.array([arange * 3, arange * 1]) # 2x40 线性正相关\n",
    "\n",
    "mean_x = np.mean(samples[0, :])\n",
    "mean_y = np.mean(samples[1, :])\n",
    "\n",
    "mean_vector = np.array([[mean_x], [mean_y]])\n",
    "scatter_matrix = np.zeros((2, 2))\n",
    "for i in range(samples.shape[1]):\n",
    "    tt = samples[:, i].reshape(2, 1) - mean_vector\n",
    "    scatter_matrix += tt.dot(tt.T)\n",
    "print('Scatter Matrix:', scatter_matrix)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Covariance Matrix \n",
    "\n",
    "协方差:\n",
    "\n",
    "$$\n",
    "Cov(X,Y)=\\frac{\\sum_{i=1}^n{(X_i-\\overline{X})(Y_i-\\overline{Y})}}{(n-1)}\n",
    "$$\n",
    "\n",
    "The maximum likelihood estimate, given n samples, for the covariance matrix of a multivariate normal distribution can be expressed as the normalized scatter matrix:\n",
    "\n",
    "$$\n",
    "C_{ML} = \\frac{1}{n} S\n",
    "$$\n",
    "\n",
    "\n",
    "协方差矩:\n",
    "\n",
    "$$\n",
    "C=\n",
    "\\begin{array}{ccc}\n",
    "cov(x,x)&cov(x,y)&cov(x,z)\\\\\n",
    "cov(y,x)&cov(y,y)&cov(y,z)\\\\\n",
    "cov(z,x)&cov(z,y)&cov(z,z)\n",
    "\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 124,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Covariance Matrix: [[1230.          410.        ]\n",
      " [ 410.          136.66666667]]\n",
      "Scatter Matrix: [[47970. 15990.]\n",
      " [15990.  5330.]]\n",
      "Unscaled covariance matrix which is same as Scatter Matrix: [[47970. 15990.]\n",
      " [15990.  5330.]]\n"
     ]
    }
   ],
   "source": [
    "print('Covariance Matrix:', np.cov(samples))\n",
    "print('Scatter Matrix:', scatter_matrix)\n",
    "print('Unscaled covariance matrix which is same as Scatter Matrix:', np.cov(samples) * 39) # n - 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Correlation Matrix\n",
    "\n",
    "相关系数:\n",
    "\n",
    "$$\n",
    "\\rho_{X,Y} = corr(X, Y) = \\frac{cov(X, Y)}{\\sigma_X \\sigma_Y}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 138,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Covariance Matrix: [[1230.          410.        ]\n",
      " [ 410.          136.66666667]]\n",
      "Correlation Matrix: [[1. 1.]\n",
      " [1. 1.]]\n",
      "Std deviation products : [[1199.25  399.75]\n",
      " [ 399.75  133.25]]\n",
      "Covariance Matrix computed from covariance : [[1.02564103 1.02564103]\n",
      " [1.02564103 1.02564103]]\n",
      "[[1.02564103 1.02564103]\n",
      " [1.02564103 1.02564103]]\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAATYAAAD+CAYAAACwcCnEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAExhJREFUeJzt3X+sX2Vhx/H3x4uFZMtG8RroSoUaGyeGWLabusVkTihQ9gdlmz9KslgXSJdFtmRGYwkJLlUS3P5gWcI2G6xUXSiOxXgXa7ryw/mHoL1mDCik9lrdaIdiWzRLkELv/eyP77nm9HLvPef2+70/+pzPqzm55/fz3Nz0k+ec55zzyDYRESV5w1JXICJi0BJsEVGcBFtEFCfBFhHFSbBFRHESbBFRnARbC5I2STokaVzS9qWuTzSTtEvSi5KeWeq6xOJLsDWQNATcC9wAXAHcLOmKpa1VtHA/sGmpKxFLI8HWbAMwbvuI7VeBPcDmJa5TNLD9LeDkUtcjlkaCrdlq4Pna8tFqXUQsUwm2iChOgq3ZMWBNbfnSal1ELFMJtmYHgHWS1kpaAWwBRpe4ThExhwRbA9ungduAfcBzwFdsH1zaWkUTSQ8AjwNvl3RU0i1LXadYPMpniyKiNGmxRURxEmwRUZwEW0QUJ8EWEcVJsM2DpG1LXYeYn/zNFkfTRwfU8/fVhySekvRbtW1bJR2upq2DqE+CbX7yn+Tck7/Z4rifuT86cAOwrpq2Af8IIOki4FPAu+m9l/0pSSv7rUyCLSL61uKjA5uBL7rnCeBCSauA64H9tk/afgnYzwC+ynJevyeYyfBFQ758zRsX4tRL6i2rz2PkXRfkwb9zSKl/sx89/xrHT06on3Nc/75f8YmTE632/d5Tpw4Cr9RW7bS9cx7FzfYxiQX5yMSCBNvla97Id/etad4xIs7Khuufb96pwYmTE3x331ta7Tu06vArtkf6LnSR5FI0oqMMTLb8NwCzfUxiQT4ykWCL6ChjXvNEq2kARoEPV72jvwP83PYL9N7Bvk7SyqrT4LpqXV8W5FI0Is4NA2qNTX104PeBYUlH6fV0vhHA9j8Be4E/AMaBl4E/rbadlPRpel/RAdhhu+8vHyfYIjrKmIkBfQTD9s0N2w18dJZtu4BdA6lIJcEW0WGTFNdhDCTYIjrLwESCLSJKkxZbRBTFwGuFfmg2wRbRUca5FI2Iwhgmysy1BFtEV/XePChTgi2is8QEfb1Hv2wl2CI6qtd5kGCLiIL0nmNLsEVEYSbTYouIkqTFFhHFMWKi0C+XJdgiOiyXohFRFCNe9dBSV2NBJNgiOqr3gG4uRSOiMOk8iIii2GLCZbbYyvytIqKVSdRqaiJpk6RDksYlbZ9h+z2Snqym70v6WW3bRG3b6CB+r7TYIjqq13nQfwRIGgLuBa6lN+DxAUmjtp/9ZVn2X9X2/wvgqtopfmF7fd8VqUmLLaKjpjoP2kwNNgDjto/YfhXYA2yeY/+bgQcG81vMLMEW0WETVqupwWqgPjT90Wrd60i6DFgLPFpbfYGkMUlPSLqpn99nSi5FIzpqnm8eDEsaqy3vtL3zLIrdAjxknzEK82W2j0l6K/CopKdt/+Aszv1LCbaIDpts3yt63PbILNuOAWtqy5dW62ayhWnji9o+Vv08Iumb9O6/9RVsuRSN6KjeS/BvaDU1OACsk7RW0gp64fW63k1JvwmsBB6vrVsp6fxqfhh4D/Ds9GPnKy22iI4y4rUBvFJl+7Sk24B9wBCwy/ZBSTuAMdtTIbcF2FONCj/lHcDnJE3Sa2jdXe9NPVsJtoiOshnYA7q29wJ7p627c9ryX89w3LeBKwdSiZoEW0RntXv49lyUYIvoKDO4Fttyk2CL6LB8aDIiimKUD01GRFl6w++VGQFl/lYR0UIGTI6Iwph5vXlwTkmwRXRYWmwRURRbabFFRFl6nQcZpSoiilLumAcJtoiO6nUe5B5bRBQmbx5ERFHy5kFEFCkjwUdEUWx4bTLBFhEF6V2KJtgiojB58yAiilLy4x6t2qGSNkk6JGlc0vaFrlRELIbepWibqfFMDRkh6SOSfirpyWq6tbZtq6TD1bR1EL9ZY4tN0hBwL3AtvRGeD0gaHcRIMhGxtAYx5sE8MuJB27dNO/Yi4FPACL1G5PeqY1/qp05tWmwbgHHbR2y/CuwBNvdTaEQsvV6v6FCrqUE/GXE9sN/2ySrM9gObzvqXqrQJttXA87Xlo9W6M0jaJmlM0thPT0xM3xwRy8zUA7ptJmB46v93NW2rnapVRgB/LOkpSQ9Jmho5vu2x8zKwzgPbO4GdACPvusANu0fEMjCPS9Hjtkf6KOrfgAdsn5L0Z8Bu4Oo+zjenNi22Y8Ca2vKl1bqIOIdN9Yq2bLHNpTEjbJ+wfapavA/47bbHno02wXYAWCdpraQV9IapH204JiLOAQPqFW3MCEmraos3As9V8/uA6yStlLQSuK5a15fGS1HbpyXdVhU2BOyyfbDfgiNiadni9ADePJgtIyTtAMZsjwJ/KelG4DRwEvhIdexJSZ+mF44AO2yf7LdOre6x2d4L7O23sIhYXgb1gO5MGWH7ztr87cDtsxy7C9g1kIpU8uZBREeV/OZBgi2iwxJsEVGUfGgyIoo0iFeqlqMEW0RH2XA6H5qMiNLkUjQiipJ7bBFRJCfYIqI06TyIiKLYuccWEcURE+kVjYjS5B5bRBQl74pGRHncu89WogRbRIelVzQiiuJ0HkREiXIpGhHFKbVXtMx2aEQ0snvB1mZqImmTpEOSxiVtn2H7xyQ9W40r+oiky2rbJiQ9WU0DGSgqLbaIDhvE4x6ShoB7gWvpDXh8QNKo7Wdru/0nMGL7ZUl/DvwN8KFq2y9sr++7IjVpsUV0mN1uarABGLd9xParwB5g85nl+DHbL1eLT9AbP3TBJNgiOsqIyck3tJqAYUljtWlb7VSrgedry0erdbO5BfhGbfmC6pxPSLppEL9bLkUjOmwenaLHbY/0W56kPwFGgPfWVl9m+5iktwKPSnra9g/6KSfBFtFVHliv6DFgTW350mrdGSRtBO4A3mv71C+rYR+rfh6R9E3gKqCvYMulaESXueU0twPAOklrJa0AtgBn9G5Kugr4HHCj7Rdr61dKOr+aHwbeA9Q7Hc5KWmwRHTaIFpvt05JuA/YBQ8Au2wcl7QDGbI8Cfwv8KvAvkgD+x/aNwDuAz0mapNfQuntab+pZSbBFdJSBycnBPKBrey+wd9q6O2vzG2c57tvAlQOpRE2CLaKrDBT65kGCLaLD8q5oRJQnwRYRZWn3Hui5KMEW0WVpsUVEUQweUK/ocpNgi+i0BFtElCaXohFRnARbRBQlD+hGRInygG5ElCe9ohFRGqXFFhFFafettXNSgi2is5TOg4goUFpsEVGcyaWuwMJIsEV0VcHPsWUwl4gOk9tNjeeRNkk6JGlc0vYZtp8v6cFq+3ckXV7bdnu1/pCk6wfxeyXYIrpsAKNUSRoC7gVuAK4AbpZ0xbTdbgFesv024B7gs9WxV9Ab1eqdwCbgH6rz9SXBFhH92gCM2z5i+1VgD7B52j6bgd3V/EPANeoNV7UZ2GP7lO0fAuPV+fqSYIvosHlcig5LGqtN22qnWQ08X1s+Wq1jpn1snwZ+Dryp5bHzls6DiK4y83ml6rjtkQWszUClxRbRZYMZCf4YsKa2fGm1bsZ9JJ0H/DpwouWx85Zgi+iwAfWKHgDWSVoraQW9zoDRafuMAlur+fcDj9p2tX5L1Wu6FlgHfLff3yuXohFdNoA3D2yflnQbsA8YAnbZPihpBzBmexT4PPAlSePASXrhR7XfV4BngdPAR21P9FunBFtElw3olSrbe4G909bdWZt/BfjALMfeBdw1mJr0JNgiOqrtw7fnogRbRJflQ5MRUZq02CKiPAm2iChK7rFFRJESbBFRGhX6ocm8eRARxUmLLaLLcikaEUVJ50FEFCnBFhHFSbBFRElEub2iCbaIrso9togoUoItIoqTYIuI0uRSNCLKU2iw5ZWqiK5yr1e0zdQPSRdJ2i/pcPVz5Qz7rJf0uKSDkp6S9KHatvsl/VDSk9W0vqnMBFtElw1m+L0m24FHbK8DHqmWp3sZ+LDtdwKbgL+TdGFt+ydsr6+mJ5sKTLBFdNiAht9rshnYXc3vBm6avoPt79s+XM3/L/Ai8OazLTDBFtFl7Vtsw5LGatO2eZRyse0XqvkfAxfPtbOkDcAK4Ae11XdVl6j3SDq/qcB0HkR01fwuM4/bHplto6SHgUtm2HTHGUXalmZvA0paBXwJ2Gp76u7e7fQCcQWwE/gksGOuyibYIjpKDO5xD9sbZy1H+omkVbZfqILrxVn2+zXg68Adtp+onXuqtXdK0heAjzfVJ5eiER22SPfYRoGt1fxW4Guvq4e0Avgq8EXbD03btqr6KXr3555pKjDBFtFli9MrejdwraTDwMZqGUkjku6r9vkg8HvAR2Z4rOOfJT0NPA0MA59pKjCXohFdtggP6No+AVwzw/ox4NZq/svAl2c5/ur5lplgi+iqfN0jIoqUYIuI0uRDkxFRnFyKRkRZBtPjuSwl2CK6LMEWESUZ5JsHy02CLaLDNFlmsiXYIroq99giokS5FI2I8iTYIqI0abFFRHkSbBFRFOeVqogoTJ5ji4gyucxkS7BFdFhabBFRloIf0G0c80DSLkkvSmocQCEizi2abDf1VYZ0kaT9kg5XP1fOst9EbbyD0dr6tZK+I2lc0oPVwC9zajOYy/30hpyPiMIsRrAB24FHbK8DHqmWZ/IL2+ur6cba+s8C99h+G/AScEtTgY3BZvtbwMnGqkfEucX0Og/aTP3ZDOyu5nfTG0KvlWrIvauBqSH5Wh0/sOH3JG2TNCZp7KcnJgZ12ohYQPMYV3R46v93NW2bRzEX1wY9/jFw8Sz7XVCd+wlJU+H1JuBntk9Xy0eB1U0FDqzzwPZOesPPM/KuCwq9JRlRmPb/U4/bHplto6SHgUtm2HTHGcXZlmbti73M9jFJbwUercYS/XnrGtakVzSiowb5gK7tjbOWI/1E0irbL1Sjur84yzmOVT+PSPomcBXwr8CFks6rWm2XAsea6pOR4CO6ykaT7aY+jQJbq/mtwNem7yBppaTzq/lh4D3As7YNPAa8f67jp2vzuMcDwOPA2yUdldTYIxER5wi3nPpzN3CtpMPAxmoZSSOS7qv2eQcwJum/6AXZ3bafrbZ9EviYpHF699w+31Rg46Wo7Zvn/WtExDlhMd48sH0CuGaG9WPArdX8t4ErZzn+CLBhPmXmHltEVxnImAcRUZwycy3BFtFleQk+IoqT4fcioiwFf90jwRbRUb0HdMtMtgRbRJdlzIOIKE1abBFRltxji4jyDOQ90GUpwRbRZbkUjYiiZMDkiChSWmwRUZwycy3BFtFlmizzWjTBFtFVJg/oRkRZhPOAbkQUqNBgy2AuEV22CAMmS7pI0n5Jh6ufK2fY532SnqxNr0yNLSrpfkk/rG1b31Rmgi2iq6busbWZ+rMdeMT2OuCRavnMqtiP2V5vez29kd9fBv69tssnprbbfrKpwARbRIdpcrLV1KfNwO5qfjdw0xz7Qm+ovW/YfvlsC0ywRXRWy8vQ/u/DXWz7hWr+x8DFDftvAR6Ytu4uSU9Jumdq/NG5pPMgoqvMfEJrWNJYbXmn7Z1TC5IeBi6Z4bg7zijStjT7SAvVSPFXAvtqq2+nF4grgJ30xhndMVdlE2wRXdb+KvO47ZHZNtreONs2ST+RtMr2C1VwvThHOR8Evmr7tdq5p1p7pyR9Afh4U2VzKRrRYbJbTX0aBbZW81uBr82x781MuwytwhBJond/7pmmAhNsEV22OPfY7gaulXQY2FgtI2lE0n1TO0m6HFgD/Me04/9Z0tPA08Aw8JmmAnMpGtFVNkws/DtVtk8A18ywfgy4tbb8I2D1DPtdPd8yE2wRXVbomwcJtoguS7BFRFEMZMyDiCiLwWV+tyjBFtFVZlE6D5ZCgi2iy3KPLSKKk2CLiLIM5OHbZSnBFtFVBjKYS0QUJy22iCjL4rxStRQSbBFdZXCeY4uI4uTNg4goTu6xRURR7PSKRkSB0mKLiLIYT0wsdSUWRIItoqvy2aKIKFKhj3tkMJeIjjLgSbea+iHpA5IOSpqUNOsQfpI2STokaVzS9tr6tZK+U61/UNKKpjITbBFd5epDk22m/jwD/BHwrdl2kDQE3AvcAFwB3CzpimrzZ4F7bL8NeAm4panABFtEh3liotXUVxn2c7YPNey2ARi3fcT2q8AeYHM1lujVwEPVfrvpjS06pwW5x/a9p04dH1o1/t8Lce4lNgwcX+pKxLyU+je7rN8T/B8v7XvYDw233P0CSWO15Z22d/Zbh5rVwPO15aPAu4E3AT+zfbq2/nVD9E23IMFm+80Lcd6lJmnM9qz3CGL5yd9sdrY3Depckh4GLplh0x225xr5fUGkVzQi+mZ7Y5+nOEZvFPgpl1brTgAXSjqvarVNrZ9T7rFFxHJwAFhX9YCuALYAo7YNPAa8v9pvK9DYAkywzc8g7ynE4sjfbIlJ+kNJR4HfBb4uaV+1/jck7QWoWmO3AfuA54Cv2D5YneKTwMckjdO75/b5xjJd6LtiEdFdabFFRHESbBFRnARbRBQnwRYRxUmwRURxEmwRUZwEW0QU5/8BazXR1QE8wwgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "print('Covariance Matrix:',np.cov(samples))\n",
    "\n",
    "std_dev_of_x1 = np.std(arange * 3)\n",
    "std_dev_of_x2 = np.std(arange * -1)\n",
    "\n",
    "std_dev_products = np.array(\n",
    "   [[std_dev_of_x1 * std_dev_of_x1, std_dev_of_x1 * std_dev_of_x2],\n",
    "    [std_dev_of_x1 * std_dev_of_x2, std_dev_of_x2 * std_dev_of_x2]]\n",
    ")\n",
    "\n",
    "pho = np.divide(np.cov(samples), std_dev_products)\n",
    "\n",
    "print('Correlation Matrix:', np.corrcoef(samples))\n",
    "print('Std deviation products :', std_dev_products)\n",
    "print('Covariance Matrix computed from covariance :', pho)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## How are these matrices related to each other?\n",
    "\n",
    "Basically, all of these matrices are calculated using the same procedure: $A^TA$. The only difference between them is how the data is scaled before the matrix multiplication is executed: \n",
    "\n",
    "- scatter matrix: no scaling ???\n",
    "\n",
    "- covariance: the mean of each variable is subtracted before multiplication\n",
    "\n",
    "- cross correlation: each variable is standardized (mean subtracted, then divided by standard deviation)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# References\n",
    "\n",
    "[Scater,Covariance,Correlation Matrix Explained](https://medium.com/@raghavan99o/scatter-matrix-covariance-and-correlation-explained-14921741ca56)\n",
    "\n",
    "[seaborn pairplot](https://seaborn.pydata.org/generated/seaborn.pairplot.html)\n",
    "\n",
    "[Visualize ML](https://machinelearningmastery.com/visualize-machine-learning-data-python-pandas)\n",
    "\n",
    "[wiki: Scatter Matrix](https://en.wikipedia.org/wiki/Scatter_matrix)\n",
    "\n",
    "[wiki: Covariance Matrix](https://en.wikipedia.org/wiki/Covariance_matrix)\n",
    "\n",
    "[wiki: Correlation Matrix](https://en.wikipedia.org/wiki/Correlation_and_dependence#Correlation_matrices)\n",
    "\n",
    "[Difference](http://www.statistics4u.com/fundstat_eng/cc_covarmat.html)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
